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The  work  under  this  contract  is  intended  to  transfer  European  data  and  thinking  about  sea-bed 
ripples  into  the  ONR  Ripples  DRI  project.  The  goal  was  to  develop  a  generic  predictor  for 
bedform  existence,  growth/decay,  height  and  spacing,  and  temporal  variability  at  a  sandy  seabed 
location  (ripples  in  currents,  waves,  and  waves-plus-currents)  as  a  function  of:  sediment 
characteristics,  water  depth,  wave  and  current  forcing,  biological  effects  and  time  history  of 
these  processes.  The  work  made  use  of  existing  data,  through  data  mining  and  interpretation,  to 
underpin  the  ongoing  SAX99  and  SAX04  collection  of  specific  sediment-acoustic  data  within 
the  ONR  program.  The  present  report  relates  to  Phase  2  of  the  project,  in  which  a  time-evolving 
predictor  for  ripple  properties  generated  by  waves  and/or  currents  (including  bio-degradation) 
was  developed.  Information  on  the  derivation  of  the  predictor  and  some  example  applications 
are  included.  The  Mark  2  Time-evolving  Ripple  Predictor  is  given  in  step-by-step  algorithmic 
form  in  Appendix  A.  The  predictor  has  been  implemented  and  tested  in  an  Excel  spreadsheet 
“Ripple  evolution  waves  +  currents  V2.0”. 
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1.  Introduction 

1.1  ONR  BAA  04-001 

This  project  has  the  title  “Development  of  a  new  marine  ripple  bedform  predictor  for 
application  in  sandy  shelf  environments”.  It  forms  part  of  the  Department  Research 
Initiative  (DRI),  ONR,  on  Critical  Benthic  Environmental  Processes  and  Modeling  at 
SAX04  (aha  Ripples  DRI).  It  refers  to  the  ONR  Long  Range  BAA  04-001,  dated  10 
Sep  2003. 

1.2  OBJECTIVES 

The  goal  was  to  develop  a  generic  predictor  for  bedform  existence,  growth/decay,  height 
and  spacing,  and  temporal  variability  at  a  sandy  seabed  location  (ripples  in  currents, 
ripples  in  waves,  ripples  in  waves  and  currents)  as  a  function  of: 

•  sediment  characteristics 

•  water  depth 

•  wave  and  current  forcing 

•  biological  effects 

•  time  history  of  the  above  processes 

The  research  made  use  of  the  knowledge  of  marine  sediment  transport  bedform 
prediction  held  by  the  Pis  and  accessed  in  the  UK/European  framework  to  the  DRI 
program.  Phase  1  of  the  project  made  use  of  existing  data,  through  data  mining  and 
interpretation,  to  underpin  the  SAX99  and  SAX04  collection  of  specific  sediment- 
acoustic  data  within  the  ONR  program.  A  Mark  1  version  of  the  ripple  predictor  was 
developed  initially  based  on  this  existing  data.  An  improved  Mark  2  predictor  was 
developed  in  Phase  2  of  the  project.  The  resulting  algorithm  was  delivered  for  use  in 
the  DRI  Ripples  program  as  Appendix  A  of  this  report. 

The  project  provides  predictive  tools  for  the  response  of  ripples  to  changes  in  wave  and 
wave-current  forcing,  including  biological  degradation. 

The  present  report  relates  to  the  Phase  2  objectives,  specifically: 

•  Visit  USA,  present  and  discuss  Mark  1  predictor,  and  requirements  for  Mark  2 

•  Test  Mark  1  ripple  predictor  with  existing  US  data  (e.g.  SAX99) 

•  Receive  data  from  SAX04  experiment 

•  Test  versus  SAX04  data,  and  adapt  predictor  accordingly 

•  Add  in  prediction  of  rates  of  growth/decay  and  migration,  and  effects  of 
bioturbation  and  “history  effect” 

•  Delivery  of  Mark  2  predictor  as  an  algorithm/subroutine 

•  Publication  in  scientific  literature  (subject  to  ONR  approval) 

The  form  of  the  ripple  predictor  should  be  geared  to  its  effect  on  sonar  performance,  and 
written  in  a  form  that  could  be  used  to  give  a  forecast  of  the  spatial  and  temporal 
variations  of  ripple  properties  over  a  sea  area,  as  a  function  of  hydrodynamic  and 
sedimentary  input  distributions. 
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There  are  a  number  of  steps  required  to  develop  a  ripple  predictor  of  the  kind  required; 

1.  Choose  or  develop  an  equilibrium  wave-generated  ripple  predictor  for  ripple  height 
and  length 

2.  Choose  or  develop  an  equilibrium  current-generated  ripple  predictor  for  ripple 
height  and  length 

3.  Devise  a  method  for  either  combining  items  1  and  2,  or  choosing  between  them, 
for  combined  wave-plus-current  conditions 

4.  Set  the  equilibrium  orientation  of  the  ripples  (i.e.  the  direction  of  the  normal  to  the 
ripple  crest-lines)  as  being  the  direction  of  the  waves  or  currents  depending  which 
is  dominant 

5.  Develop  a  method  of  predicting  the  rate-of-change  of  wave-ripple  geometry  as  a 
function  of  the  wave  and  sediment  characteristics 

6.  Develop  an  analogous  method  for  the  rate-of-change  of  current-ripple  geometry 

7.  Develop  a  mathematical  method  for  including  biological  effects  on  ripple 
geometry 

8.  Combine  items  1  to  7  into  a  single  algorithm 

9.  Test  the  algorithm  against  lab  and  field  data. 

The  above  steps  are  described  in  detail  in  Sections  3  to  7  of  this  report. 
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2.  Behaviour  of  sea-bed  ripples 

2.1  RIPPLES  AND  SONAR  PROPAGATION 

The  behaviour  of  ripples  (or  other  types  of  bedform),  as  they  might  affect  sonar 
propagation,  can  be  considered  in  three  stages: 

1.  Prediction  of  the  existence  of  different  types  of  bedform  under  different 
hydrodynamic  conditions  (waves,  currents,  water  depths)  and  sedimentary 
conditions  (characteristic  grain-sizes,  density,  mineralogy). 

2.  Prediction  of  the  characteristic  descriptors  of  ripples  in  equilibrium  with  the 
hydrodynamic  conditions  (i.e,  assuming  that  the  conditions  have  lasted  sufficiently 
long  and  been  sufficiently  constant  that  the  ripples  have  evolved  fully). 

3.  Prediction  of  the  characteristic  descriptors  of  ripples  in  conditions  that  are  not  in 
equilibrium ,  and  how  fast  they  respond  (tidally  varying  currents,  wind-driven 
varying  currents,  varying  waves  in  storms,  biological  effects). 

Phase  1  of  the  project  delivered  items  1  and  2,  and  Phase  2  delivered  item  3. 

The  term  bedform  is  used  here  to  encompass  any  kind  of  deviation  from  a  flat  bed,  of 
which  ripples  are  the  most  prevalent,  but  can  also  include  dunes,  mega-ripples, 
hummocks  and  sandwaves.  We  distinguish  ripples  as  small-scale  bedforms,  having 
spacing  of  order  lm  and  heights  of  order  0.1m,  which  enhance  the  penetration  of  sonar 
signals  into  the  subsurface.  Ripples  (and  other  bedforms)  may  be  two-dimensional 
(2D),  in  which  case  the  ripple  crests  form  straight  or  gently  waving  lines  that  are  very 
much  longer  than  the  wavelength  (i.e.  spacing  perpendicular  to  the  crest-line),  or  three- 
dimensional  (3D),  in  which  case  the  shape  of  an  individual  ripple  can  be  traced  (along 
the  crest-line)  for  only  a  short  distance.  Both  2D  and  3D  ripples  (and  bedforms)  are 
encountered  on  the  sea-bed,  but  their  effects  on  sonar  propagation  will  be  different. 

We  consider  the  following  properties  of  ripples  that  might  affect  sonar  propagation: 

•  height  (from  trough  to  crest) 

•  wavelength  (i.e.  spacing  perpendicular  to  the  crest-lines) 

•  crest-length  (i.e.  spacing  along  the  crest-lines) 

•  orientation  (with  respect  to  the  North,  or  sonar  source) 

•  slope  with  respect  to  horizontal  (as  seen  by  sonar  source) 

•  shape 

•  grain-size. 

For  each  of  these,  the  effect  we  might  expect  ripples  to  have  on  sonar  are  as  follows: 

•  height:  the  sonar  reflection  or  transmission  will  increase  with  ripple  height 

•  spacing:  the  interaction  between  ripples  and  sonar  will  depend  on  the  relation 
between  the  sonar  wavelength  (projected  onto  the  bed)  and  the  spacing  of  the 
ripples,  being  strongest  if  these  are  matched.  This  interaction  will  be  stronger  if 
the  periodicity  of  the  ripples  is  “sharp”  rather  than  diffuse 

•  crest-length:  long-crested  (2D)  ripples  will  have  little  impact  on  a  sonar  beam 
directed  along  the  crest-line,  whereas  short-crested  (3D)  ripples  will  impact  a  sonar 
beam  even  if  it  is  directed  along  the  (less  well-defined)  crest-line 
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•  orientation:  the  strength  of  the  interaction  will  be  greatest  if  the  sonar  source  faces 
the  steepest  slopes  of  the  ripples  directly  (in  azimuth),  and  will  be  progressively 
weaker  as  the  ripple  crest-lines  form  greater  azimuths  with  the  sonar  direction 

•  slope:  the  reflection/transmission  of  sonar  will  be  greatest  if  the  face  of  the  ripple- 
slopes  is  most  nearly  perpendicular  to  the  sonar  beam.  Hence  asymmetric  current- 
induced  ripples  will  have  greatest  interaction  with  a  sonar  beam  directed  towards 
the  steeper  (lee)  slope,  and  a  weaker  interaction  if  it  is  directed  towards  the  gentler 
(stoss)  slope 

•  shape:  reflection  or  transmission  will  be  greater  for  sharp-crested  (triangular) 
cross-section  ripples  than  for  round-crested  ripples 

•  grain-size:  the  absorption  or  reflection  of  sonar  will  also  depend  on  the  grain-size 
of  the  sediment  forming  the  ripple. 

In  this  project  the  main  attention  was  focussed  on  the  most  important  geometric  ripple 
properties:  height,  wavelength  and  orientation.  We  use  the  words  ripple  length, 
wavelength,  and  spacing  interchangeably. 

2.2  RESULTS  FROM  PHASE  1 

The  following  tasks  were  achieved  in  Phase  1  of  the  Research. 

HR  Wallingford: 

•  held  discussions  with  European  (and  other)  ripples  researchers 

•  compiled  a  digest  of  European  source  publications 

•  assembled  a  data-base  of  ripple  geometries  (height,  wavelength)  and  driving 
conditions  (waves,  currents,  sediment  size) 

•  performed  statistical  analyses  of  ripple  heights  and  wavelengths 

•  produced  a  graphical  display  in  Excel  of  colour-coded  bed  elevations  and  cross- 
sections  in  two  horizontal  dimensions 

•  presented  the  results  in  HR  Wallingford  Report  TR150  (Soulsby  and  Whitehouse, 
2005). 

Further  details  of  all  these  results  were  presented  in  HR  Wallingford  Report  TR150,  and 
are  not  repeated  here. 

2.3  RESULTS  FROM  PHASE  2 

At  the  end  of  Phase  2,  the  following  tasks  have  been  achieved. 

HR  Wallingford: 

•  devised  new  formulas  to  predict  equilibrium  wave-generated  ripple  heights  and 
lengths  as  functions  of  wave  and  sediment  properties,  and  tested  them  against  the 
data-base 

•  presented  results  at  a  SAX04  workshop  in  Seattle  (May  2005) 

•  submitted  an  abstract  (now  accepted)  to  ICCE  conference,  San  Diego  Sept  2006, 
on  wave-generated  ripple  predictor  (copied  to  ONR) 

•  devised  alternative  new  formulas  for  predicting  equilibrium  wave-generated 
ripples,  and  used  Phase  1  data-base  to  choose  the  most  successful 

•  devised  new  formulas  for  current-generated  ripples  as  functions  of  current  speed 
and  sediment  properties,  and  tested  against  data-base 
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•  brought  together  these  formulas  for  equilibrium  ripple  heights  and  wavelengths 
(and  orientations)  to  form  the  “Mark  1”  ripple  predictor 

•  assembled  information  about  existing  intercomparisons  of  ripple  predictors  (wave, 
current  and  combined)  against  existing  lab  and  field  data-sets,  and  summarised  this 
in  Excel  spreadsheets 

•  produced  a  pro-forma  spreadsheet  for  collecting  SAX04  data  from  experimenters, 
and  submitted  it  to  ONR 

•  devised  a  new  method  for  predicting  the  time-development  of  ripple  heights, 
wavelengths  and  orientations 

•  encapsulated  the  equilibrium  predictors  in  the  time-development  model  to  form  the 
“Mark  2”  ripple  predictor,  as  algorithms  and  an  Excel  spreadsheet 

•  added  in  a  term  to  simulate  bio-degradation  of  ripples 

•  tested  the  time-development  model  for  wave-generated  ripple  evolution  against  an 
existing  time-series  of  wave  data  measured  by  University  of  Aberdeen  (UK)  in  a 
laboratory  oscillating  water  tunnel  (effectively  at  full  scale),  and  obtained  good 
agreement 

•  tested  the  time-development  model  for  current-generated  ripple  evolution  against 
an  existing  time-series  of  ripple  observations  measured  by  HR  Wallingford  for 
tidally  oscillating  current  speed  in  a  laboratory  flume  (effectively  at  full  scale),  and 
obtained  good  agreement 

•  tested  the  applicability  of  the  time-development  model  for  wave-plus-current¬ 
generated  ripple  evolution  against  data  measured  by  HR  Wallingford  off  the  UK 
coast  (but  without  ripple  observations). 

At  the  end  of  the  present  contract  we  have  developed  and  tested  a  Mark  2  predictor  for 
the  time-evolution  of  ripple  heights,  wavelengths  and  orientations,  driven  by  any  time- 
series  of  field  inputs  of  waves  and  currents  (e.g.  hourly  values  for  several  weeks  or 
months)  for  sand-sized  sediments.  This  predictor  has  the  capability  to  not  only  make 
reasonably  accurate  predictions  of  ripple  height  and  wavelength  under  steady 
wave/current  conditions,  but  also  allows  the  ripples  to  grow  or  decay  at  a  rate  controlled 
by  the  wave  and  current  inputs,  and  with  the  option  of  including  biological  degradation 
of  ripples.  Thus  for  example,  if  the  wave/current  conditions  decrease  below  the 
threshold  of  motion,  the  predicted  ripples  become  “frozen”  (sometimes  called  relict 
ripples,  or  hysteresis),  then  re-adjust  their  height,  wavelength  and  orientation  once 
conditions  again  become  intense  enough.  Gradual  ripple  washout  by  very  intense 
conditions,  or  decay  due  to  biological  processes,  are  included,  with  subsequent  re¬ 
growth  when  suitable  conditions  return.  We  believe  this  is  the  first  ripple  predictor  to 
have  all  these  capabilities. 

However,  we  have  not  been  able  to  test  the  model  against  SAX99  data  as  intended, 
because  it  proved  to  be  unsuitable  (only  “frozen”  ripples  were  observed,  and  the 
wave/current  conditions  forming  them  had  not  been  measured).  We  have  also  not  been 
able  to  test  the  model  against  the  SAX04  data-set  as  intended,  because  the 
experimenters  wish  to  distil  and  publish  their  own  results  from  it  before  releasing  it  for 
third  party  use. 
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3.  Existing  ripple  geometry  predictors 

3.1  EXISTING  RIPPLE  PREDICTORS 

Before  devising  new  ripple  predictors,  a  study  was  made  to  see  if  there  were  existing 
predictors  that  met  all  the  present  requirements.  This  was  based  on  examination  of 
intercomparison  studies  reported  in  the  scientific  literature.  The  statistics  of  the  study 
are  tabulated  in  Tables  1  -  5  in  terms  of  ripple  predictors  and  the  ripple  data-sets  they 
were  intercompared  with.  It  should  be  noted  that  the  predictors  and  data-sets  listed  in 
these  tables  are  only  those  used  in  this  set  of  intercomparison  studies.  They  are  not 
exhaustive,  and  we  have  made  use  of  other  data-sets  in  both  HR  report  TR  150  (Soulsby 
and  Whitehouse,  2005)  and  the  present  report. 

The  information  contained  in  Tables  1  -  5  is  as  complete  and  accurate  as  we  can  make  it 
without  having  read  every  reference.  Some  of  the  collected  information  on  data-sets 
and  predictors  has  been  extracted  from  third  party  sources  and  we  have  relied  on  the 
information  presented.  Overall,  the  tables  are  intended  to  provide  an  indication  of  the 
situation  at  the  present  time  based  on  our  experience  of  the  literature. 

A  list  of  existing  ripple  predictors  has  been  compiled  from  the  study  of  the  literature. 
These  are  variously  for  application  to  height  and/or  wavelength  of  ripples,  generated  by 
waves,  (W),  currents  (C)  or  both  together  (W+C),  and  for  either  equilibrium  (E)  or 
transient  (T)  conditions.  The  methods  (and  intercomparison  papers),  dating  from  1980 
up  to  2005,  are  listed  in  Table  1.  A  code  is  assigned  to  each  method  (e.g.  GM82  for 
Grant  and  Madsen,  1982)  for  use  in  subsequent  tables  and  text.  A  total  of  39  methods 
and  intercomparisons  are  listed,  with  28  for  waves  alone,  5  for  currents  alone,  and  6  for 
combined  waves  and  currents.  Of  the  predictors,  18  are  for  waves  alone,  5  for  currents 
alone,  and  5  for  combined  waves  and  currents.  Only  4  of  the  methods  treat  the  transient 
(time-varying)  case,  the  remainder  being  for  equilibrium  (steady)  conditions.  A  more 
detailed  analysis  of  9  of  the  wave-only  predictors  is  given  in  Section  3.4. 

3.2  EXISTING  DATA-SETS  USED  TO  TEST  PREDICTORS 

A  list  of  data-sets  has  been  compiled  based  on  those  used  to  test  the  predictors  listed  in 
Table  1.  This  is  not  intended  to  be  an  exhaustive  catalogue  of  ripple  data-sets,  but  only 
those  used  for  testing  and  intercomparison  of  predictors.  They  are  categorised  as 
generated  by  waves  (W),  currents  (C)  or  both  together  (W+C);  for  equilibrium  (E)  or 
transient  (T)  conditions;  and  from  the  laboratory  (L)  or  the  field  (F).  The  data-sets, 
dating  from  1939  up  to  2005,  are  listed  in  Table  2.  As  with  the  predictors,  codes  are 
assigned  for  subsequent  reference.  A  total  of  83  data-sets  is  listed,  of  which  50  are  for 
waves  alone,  22  are  for  currents  alone,  and  1 1  are  for  combined  waves  and  currents.  75 
data-sets  were  for  equilibrium  conditions  and  8  for  transient  conditions.  71  data-sets 
were  from  the  lab,  and  14  from  the  field. 

3.3  INTERCOMPARISON  OF  INTERCOMPARISONS 

A  total  of  24  intercomparison  papers  or  books  have  been  identified.  These  are  shown  in 
Table  3,  using  the  codes  assigned  in  Table  1.  The  predictors  tested  are  shown  with  an  x 
against  each  paper.  The  most  widely  tested  predictors  are:  Nielsen  (1981)  [13  times], 
Wiberg  and  Harris  (1994)  [9  times],  Grant  and  Madsen  (1982)  [7  times]  and  Mogridge 
et  al  (1994)  [7  times],  with  all  the  others  being  tested  3  or  less  times.  It  is  noteworthy 
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that  the  wave-only  predictors  have  been  subjected  to  far  more  testing  than  current-only 
or  W+C  predictors. 

The  most  comprehensive  intercomparisons  in  terms  of  number  of  predictors  tested  are 
those  of  Traykovski  et  al  (1999),  Foti  and  Faraci  (2003),  and  Williams  et  al  (2004),  each 
testing  6  predictors.  All  other  intercomparisons  were  of  5  or  less  predictors. 

The  data-sets  against  which  intercomparisons  have  been  made  are  listed  in  Table  4, 
using  the  codes  assigned  in  Tables  1  and  2.  The  most  widely  used  data-sets  are:  Inman 
(1957)  [10  times],  Carstens  et  al  (1969)  [9  times],  Kennedy  and  Falcon  (1965)  [7  times], 
Mogridge  and  Kamphuis  (1972)  [6  times],  and  Dingier  (1974)  [5  times],  with  all  the 
others  being  tested  4  or  less  times. 

The  most  comprehensive  intercomparisons  in  terms  of  number  of  data-sets  used  are 
those  of  Baas  (1993)  [17  data-sets],  Mogridge  et  al  (1993)  [16  data-sets],  ABP  (2004) 
[13  data-sets],  and  Mogridge  et  al  (1994)  [13  data-sets].  All  the  other  intercomparisons 
used  9  or  less  data-sets. 

The  correspondence  between  predictors  and  data-sets  in  intercomparisons  is  shown  in 
Table  5,  where  all  the  intercomparisons  are  taken  together.  Overall,  it  can  be  seen  that 
many  of  the  predictors  have  been  tested  against  a  wide  range  of  data-sets  (taking  all 
intercomparisons  into  account),  although  some  have  only  been  tested  against  the 
originator’s  own  data.  The  predictors  which  have  been  compared  with  the  largest 
number  of  data-sets  overall  are  Nielsen  (1981)  [53  data-sets],  Mogridge  et  al  (1994)  [36 
data-sets]  and  Wiberg  and  Harris  (1994)  [24  data-sets].  All  the  others  have  been  tested 
against  18  or  less  data-sets. 

No  attempt  has  been  made  here  to  rank  the  predictors,  partly  because  intercomparisons 
rarely  announce  an  unequivocal  “winner”  (unless  it  is  their  own  method!),  and  partly 
because  the  older  methods  have  been  more  frequently  tested,  whereas  one  would  expect 
more  recent  methods  to  be  better  due  to  improvements  in  data  and  understanding. 

However,  individual  intercomparisons  have  been  of  a  relatively  small  sub-set  of  the 
available  predictors,  tested  against  a  relatively  small  sub-set  of  the  available  data-sets. 
Furthermore,  the  selected  sub-sets  of  data  vary  widely  between  the  intercomparisons.  A 
strong  case  could  be  made  for  a  concerted  effort  by  a  group  of  researchers  to  assemble 
all  the  (easily  available)  data-sets,  and  perform  a  consistent  and  unbiased  test  of  the 
available  ripple  predictors. 

3.4  ANALYSIS  OF  EXISTING  WAVE-RIPPLE  PREDICTORS 

As  a  lead-in  to  choosing  or  devising  a  new  wave-generated  ripple  geometry  predictor, 
we  make  a  comparative  review  of  nine  existing  methods  for  predicting  the  height  and 
wavelength  of  wave-induced  ripples  which  are  summarised  in  Table  6.  They  are 
referenced  and  compared  in  the  ten  books  and  intercomparison  papers  listed  in  Table  7. 
Since  most  of  the  methods  are  given  in  full  in  either  Foti  and  Faraci  (2003)  [6  methods] 
or  Grasmeijer  and  Kleinhans  (2004)  [3  methods],  the  full  formulae  are  not  repeated 
here.  Foti  and  Faraci  (2003)  additionally  give  Matlab  codes  to  compute  the  6  methods 
they  reviewed. 

The  various  methods  were  originally  written  using  a  wide  variety  of  non-dimensional 
groupings  of  the  input  parameters.  However,  they  can  be  re-cast  mathematically  to  a 
smaller  number,  which  are  listed  below. 
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Output  parameters 
r|  =  ripple  height 

X  =  ripple  wavelength 

Input  parameters 

d50  =  median  grain  diameter  of  sediment 

ps  =  density  of  sediment 

p  -  density  of  water 

v  =  kinematic  viscosity  of  water 

g  =  acceleration  due  to  gravity 

Uw  =  amplitude  of  near-bed  wave-induced  orbital  velocity 
T  =  wave  period. 

Uw  and  T  are  representative  of  regular,  sinusoidal  waves.  O’Donoghue  et  al  (2005) 
showed  from  comparative  laboratory  tests  that  in  regular,  asymmetric  waves  the  best 
representation  for  Uw  is  given  by  Umax,  the  maximum  (usually  onshore-directed)  velocity 
(under  the  wave  crest),  and  for  irregular,  asymmetric  or  symmetric  waves  the  best 
representation  for  Uw  is  given  by  Ui/io,  the  mean  of  the  highest  one-tenth  velocities.  In 
irregular  waves,  the  peak-period  Tp  gives  the  best  representation  of  T. 

Derived  parameters 


A 

=  UwT/(2tc) 

=  amplitude  of  near-bed  wave  excursion 

(3.1) 

A 

“  A/d50 

(3.2) 

T 

u2 

H' 

=  wave  mobility  parameter 

(3.3) 

9'w 

\fu 2 

2  J  w  w 

-  skin-friction  Shields  parameter 

(3.4) 

g{s-l)dso 

fw 

=  grain-related  wave  friction  factor  (function  of  A  and  Rew) 

(3.5) 

Rew 

=  ■ 

V 

—  wave  Reynolds  number 

(3.6) 

s 

=  Ps/p 

=  density  ratio 

(3.7) 

D. 

= 

sfc  -  0 

V2 

X 

^50 

=  dimensionless  grain-size 

(3.8) 

X 

^50 

=  period  parameter 

(3.9) 

g{s-  \y 

1 

The  nine  methods  can  be  re-cast  to  give  non-dimensional  expressions  for  the 
dimensionless  ripple  height  r|/A  and  wavelength  A/A.  The  non-dimensional 
combinations  of  input  parameters  for  each  method  are  given  in  Table  6.  They  fall  into 
three  families: 

a.  Functions  of  A.  WH94  can  be  re-cast  to  give  r\/A  and  A/A  as  functions  of  only  A 

(Malarkey  and  Davies,  2003).  T99  also  gives  functions  of  only  A.  They  can  be 

regarded  as  “kinematic”  methods,  since  A  is  a  ratio  of  lengths  (or  velocities),  but  g  and  s 
are  not  included.  M94  gives  functions  containing  %  as  well  as  A,  introducing  a 
dependence  on  wave  period  combined  with  g  and  s. 

b.  Functions  of  T.  N81,  VR89  and  GK04  give  r|/A  and  A/A  as  functions  of  only  'F. 
FF02  gives  functions  containing  A  and  Rew  as  well  as  'F,  introducing  a  dependence  on 
viscosity.  They  can  be  regarded  as  “semi-dynamic”  methods,  because  they  include  g 


TR  154 


9 

UNCLASSIFIED 


R.  2.0 


Prediction  of  Ripple  Properties  in  Shelf  Seas 
Mark  2  Predictor  for  Time  Evolution 


UNCLASSIFIED 


Wallingford 


and  s,  but  do  not  include  fw  which  is  necessary  to  describe  the  tractive  force  on  sediment 
grains. 

c.  Functions  of  0'w.  GM82  and  Ma93  give  r|/A  and  A/A  as  functions  of  0'w  and  D*. 
They  can  be  regarded  as  “dynamic”  methods,  because  the  Shields  parameter  0'w 
represents  the  ratio  of  forces  acting  on  sediment  grains.  D*  introduces  a  dependence  on 
viscosity. 

Regarding  the  applicability  of  these  input  parameters,  the  following  observations  can  be 
made  about  wave-induced  ripples. 

1.  For  very  small  wave-induced  velocities  (Uw),  the  sediment  is  below  the  threshold 
of  motion.  Ripples  do  not  form  on  an  initially  flat  bed.  (In  the  sea,  “frozen”  or 
“relict”  ripples  remain  from  previous  more  energetic  conditions,  unless  biological 
activity  has  flattened  them.) 

2.  For  intermediate  wave-induced  velocities,  ripples  form.  At  larger  velocities  the 
ripple  heights  reduce. 

3.  For  large  wave-induced  velocities,  ripples  are  washed  out  leaving  a  flat  bed  with 
oscillatory  sheet  flow. 

4.  For  grainsizes  less  than  about  0.06mm,  ripples  are  low  and  indistinct. 

5.  For  grainsizes  between  0.06  and  about  0.7-1. 0mm,  ripples  are  well-defined. 

6.  For  grainsizes  larger  than  about  0.7- 1.0mm  ripples  are  either  indistinct  or  do  not 
form. 

Since  observations  1  to  3  relate  to  velocities,  and  4  to  6  relate  to  grainsizes,  they  cannot 
be  satisfied  only  by  a  single  parameter  such  as  A,  *¥  or  9'w,  all  of  which  contain  both 
velocity  and  grainsize.  Hence  it  would  seem  to  be  necessary  for  q/A  and  A/A  to  be 
formulated  as  independent  functions  of  parameters  such  as  A,  or  0'w  to  satisfy 
observations  1  to  3,  and  D*  or  %  to  satisfy  observations  4  to  6.  The  methods  of  GM82 
(containing  0'w  and  D*)  and  M94  (containing  A  and  %)  seem  to  satisfy  this  requirement 
best.  The  method  of  Ma93  superficially  appears  to  meet  the  requirements  but  in  fact  is  a 

function  of  only  the  single  parameter  (o'v  /D^),  so  the  dependencies  on  velocity  and 
grainsize  are  linked. 

All  the  methods  are  essentially  empirical  curves  fitted  to  lab  and  field  data.  Hence  it 
might  be  expected  that,  despite  their  dependence  on  different  parameters,  they  would 
give  similar  predictions  for  q  and  X  for  given  input  values.  We  have  tested  this  by 
plotting  all  9  methods  for  a  wide  range  of  inputs.  We  made  use  of  the  Matlab  codes 
provided  by  Foti  and  Faraci  (2003)  for  6  of  the  methods,  but  with  a  spurious  factor  of 
SQRT(2)  removed  from  the  orbital  velocities.  The  outputs  were  checked  for  the  method 
WH94  against  curves  in  the  original  paper,  but  the  other  methods  have  not  been 
checked.  The  remaining  3  methods  were  coded  in  Matlab  by  editing  the  codes  given  by 
Foti  and  Faraci  (2003). 

A  randomly  distributed  set  of  1,000  inputs  was  generated  by  fixing  the  values  p  = 
1027kg  m'3,  v  =  1.36  x  10'6mV*  and  ps  =  2650kg  m  3  (quartz  sand  in  sea  water  at  10°C, 
35ppt),  and  randomly  selecting  values  of  the  other  input  parameters  within  the  ranges 
0.2<Uw<2.0m.s'1,  4<T<16s,  and  0.06<d5o<2.0mm.  The  1,000  inputs  thus  cover  most  of 
the  conditions  likely  to  be  encountered  for  sand  in  the  sea.  The  same  1,000  inputs  were 
fed  into  the  9  prediction  methods  listed  in  Table  6,  and  the  output  values  of  q/A  and  A/A 
plotted  against  each  of  the  main  independent  variables  A,  and  0'w.  Figures  1, 2  and  3 
show  that  the  9  methods  differ  very  widely  in  their  predictions.  The  methods  that 
depend  on  only  A  plot  as  single  curves  when  A  is  the  abscissa  (Figure  1),  as  would  be 
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expected,  whereas  other  methods  scatter  widely.  The  same  is  true  of  the  plots  versus  T 
(Figure  2)  and  0\v  (Figure  3).  Irrespective  of  whether  r\  or  X  is  plotted,  and  whether  A, 
T  or  0'w  is  the  abscissa,  the  methods  vary  by  factors  of  between  10  and  1000  in  their 
predictions.  Thus  it  is  not  the  case  that  the  methods  give  similar  predictions. 

Although  there  is  already  a  wide  selection  of  ripple  predictors  available,  it  is  clear  that 
they  have  a  large  measure  of  disagreement.  Furthermore,  while  several 
intercomparisons  have  been  made  testing  predictors  against  data,  there  is  no  consistency 
in  either  the  predictors  tested  or  the  data-sets  used.  There  is  thus  no  predictor  that  is 
unequivocally  superior  to  the  others.  For  this  reason,  (and  despite  the  large  number 
already  proposed)  we  have  developed  our  own  wave-generated  ripple  predictor,  based 
partly  on  the  data-set  assembled  for  this  project,  and  partly  on  experiences  with  existing 
predictors. 
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4.  New  wave-only  ripple  geometry  predictor 

4.1  NEW  WAVE-RIPPLE  EQUILIBRIUM  PREDICTOR 

The  first  question  to  address  when  devising  a  new  predictor  for  wave-generated  ripple 
geometry  is  which  out  of  A,  \| /  and  0'w  (or  some  combination  of  them)  is  the  best  choice 
of  independent  variable.  This  has  been  tested  using  the  data-base  assembled  in  Phase  1 
to  see  which  parameterisation  makes  the  data  collapse  most  closely  onto  single  curves. 
Figures  4a,  b  and  c  show  plots  of  r\/A  versus  A,  vj/  and  0'w  respectively,  and  Figures  5a,  b 
and  c  show  similar  plots  for  A/A.  All  the  plots  show  considerable  scatter,  but  the  plots 
against  A  appear  to  cluster  the  data  more  tightly  than  those  against  \|/,  with  the  plots 
against  0'w  giving  the  greatest  scatter.  On  these  grounds,  we  choose  to  formulate  the 
predictors  for  equilibrium  wave-generated  ripple  height  and  wavelength  as  functions  of 
the  excursion  parameter  A. 

This  is  the  same  choice  as  that  used  by  Wiberg  and  Harris  (1994),  (re-stated  by 
Malarkey  and  Davies  (2003)  in  non-iterative  form).  An  alternative  case  could  be  made 
for  expressing  the  wavelength  in  terms  of  the  kinematic  ratio,  A,  but  with  the  ripple 
steepness  as  a  function  of  the  dynamical  ratio  0\v.  Dependence  of  steepness  r\/X  as  a 
function  of  0'w  was  advocated  by  Grant  and  Madsen  (1982)  (with  additional  dependence 
on  D*),  and  by  Nielsen  (1992).  However,  on  the  basis  of  the  data  collapse  shown  in 
Figures  4a  and  5a,  we  will  adhere  to  the  choice  of  A  for  the  functional  dependence  of 
both  r|  and  A,. 

We  will  follow  the  general  observations  made  by  Wiberg  and  Hams  (1994),  who  in  tum 
followed  Clifton  (1976): 

•  For  values  of  A  less  than  about  750,  A,  is  proportioned  to  A  (orbital  ripples),  with 
A/A  approximately  unity 

•  For  intermediate  values  of  A,  A/d50  (not  plotted)  reaches  a  peak  value  of  rather 
more  than  1000,  and  then  declines  with  increasing  A  (sub-orbital  ripples) 

•  For  values  of  A  larger  than  about  3000,  A,  is  proportional  to  d50  (anorbital  ripples) 
with  XJdso  being  around  500 

•  For  values  of  A  less  than  about  1500,  the  ripple  steepness  r\/\  is  roughly  constant, 
lying  between  0. 1  and  0.2 

•  As  A  increases  above  1500,  t^/A,  decreases  rapidly  with  A  as  the  ripples  are  washed 
out 

•  For  values  of  A  greater  than  about  6000,  ripples  cannot  exist  and  the  bed  is  flat  (or 
gently  undulating)  with  sheet  flow  of  sediment. 

The  following  expressions  (which  replace  those  presented  at  the  SAX04  workshop  in 
Seattle  in  May  2005)  for  ripple  wavelength  and  steepness  follow  the  general  principles 
outlined  above,  with  values  of  the  coefficients  based  on  fits  to  the  data-base  and 
similarities  with  WH94: 


X 


A 


1  +  1.87x1 0  .All  -  expb-  (2.0  x  1 0 


'4f 


(4.1) 


n 


=  0.15 


1 


exp  {-(5000  /A)35}] 


(4.2) 
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It  was  found  that  alternative  expressions  giving  q/A  directly  as  a  function  of  A  yielded 
curves  of  steepness  r\/X  versus  A  that  had  unrealistic  behaviour.  It  was  therefore  decided 
(in  common  with  many  other  methods)  to  parameterise  r\/X  in  terms  of  A,  and  then 
calculate  ripple  height  q  from: 


n=A.n 

A  A  X 


(4.3) 


An  additional  constraint  is  that  ripples  can  only  form  and  evolve  if  the  wave-induced 
stress  exceeds  the  threshold  of  motion  of  the  sediment.  This  is  expressed  by  comparing 
the  actual  Shields  parameter  0'w  with  the  threshold  value  0cr.  Then: 

•  If  0  w  ^  0cr,  ripples  take  pre  existing  values  of  t|  and  X 

•  If  0'w  >  0cr,  equilibrium  values  of  q  and  X  are  given  by  Eqs  (4.1),  (4.2)  and  (4.3). 

Tests  of  Eqs  (4.2,  3)  and  (4.1)  against  lab  and  field  data  from  the  data-base  are  shown  in 
Figures  6  and  7  respectively.  In  both  cases,  despite  a  large  amount  of  scatter  in  the  data 
points,  the  prediction  curves  pass  through  the  centre  of  the  data  reasonably  well.  It 
should  be  borne  in  mind  that  although  field  data  might  be  considered  more  relevant  than 
lab  data,  there  is  always  uncertainty  about  whether  ripples  in  the  sea  are  truly  in 
equilibrium  with  the  instantaneous  wave  and  current  conditions,  and  also  that  they  are 
rarely  due  entirely  to  waves  alone. 


4.2  TIME-EVOLUTION  OF  WAVE-RIPPLES 

As  pointed  out  above,  since  ripples  on  the  sea-bed  are  constantly  varying  in  response  to 
the  varying  wave  conditions,  the  prediction  of  the  rate  of  response  of  the  ripples  is  just 
as  important  as  prediction  of  the  equilibrium  geometry.  Effectively,  the  ripples  are 
constantly  trying  to  catch  up  with  the  driving  conditions.  In  addition,  if  the  wave- 
induced  orbital  velocities  fall  below  the  threshold  of  sediment  the  pre-existing  ripples 
remain  “frozen”.  The  time-evolution  of  ripple  height,  wavelength  and  orientation  all 
need  to  be  considered. 

Experiments  on  the  evolution  of  wave-induced  ripples  have  been  performed  by 
Doucette  and  O’Donoghue  (2005a,  2005b)  in  a  large  oscillating  water  tunnel  at  the 
University  of  Aberdeen,  UK.  Smaller  scale  experiments  were  also  undertaken  by  Smith 
and  Sleath  (2005)  in  an  oscillating  tray  rig  at  the  University  of  Cambridge,  UK.  For  the 
present  purposes,  we  concentrate  on  the  large-scale  experiments  at  Aberdeen,  referred  to 
as  DO05.  These  were  performed  using  sand  with  d50  =  0.44mm,  and  with  velocities 
simulating  the  full-scale  near-bed  flow  produced  by  asymmetric  regular  and  irregular 
waves. 

DO05  presented  plots  showing  how  ripple  height  and  wavelength  evolve  with  time, 
when  steady  wave  conditions  are  applied  to  either  an  initially  flat  bed  (ripple  growth),  or 
a  rippled  bed  produced  earlier  by  different  wave  conditions  (ripple  evolution).  The 
evolution  of  ripples  from  both  initially  larger,  and  initially  smaller,  ripples  were 
measured.  DO05  interpreted  the  results  by  fitting  an  empirical  relationship  between  a 
rate-of-change  parameter  and  the  wave  mobility  parameter  \j/.  Although  their  proposed 
formulation  gives  a  good  fit  to  their  observed  rates  of  change  of  ripple  height,  it  is 
difficult  to  apply  for  general-purpose  use  in  field  conditions,  as  its  extrapolated 
behaviour  in  unrealistic.  In  particular,  extrapolation  to  \\f  =  0  (no  waves)  would  still 
predict  some  ripple  evolution,  and  extrapolation  to  large  \\ /  predicts  excessively  rapid 
evolution.  Consequently,  a  more  suitable  formulation  has  been  developed  by  the 
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present  first  author  (RLS)  together  with  Dr  O’Donoghue  and  Dr  Doucette,  as  follows. 
DO05  (and  Smith  and  Sleath,  2005)  showed  that  the  evolution  of  ripple  height  r|(t)  from 
an  initial  value  of  r|0  to  a  new  equilibrium  value  r|eq  corresponds  closely  to  an 
exponential  relaxation  given  by: 


(4.4) 


where  T  is  wave  period,  and  p  is  a  coefficient  that  governs  the  rate  of  change  (small  for 
small  waves,  large  for  large  waves).  Eq  (4.4)  applies  for  a  step  change  from  one  steady 
wave  condition  (or  flat  bed)  to  another  steady  condition.  A  more  general  expression, 
better  suited  to  application  to  constantly  varying  field  conditions,  is  given  by: 


(4.5) 


Eq  (4.4)  is  a  solution  of  Eq  (4.5)  for  the  special  case  where  q  =  rjo  at  t  =  0,  and  p,  T  and 
r|Cq  are  independent  of  t.  The  method  of  solving  the  discretised  form  of  Eq  (4.5)  is 
described  in  Section  6.  The  following  expression  for  P  gives  a  good  fit  to  the  DO05 
data,  and  has  the  plausible  asymmetry  behaviour  of  tending  to  zero  for  zero  waves,  and 
tending  to  a  (large)  constant  for  very  large  waves: 


|3  =  2.996\|/i07/(21700  +  v|/107) 


(4.6) 


The  application  of  Eqs  (4.5)  and  (4.6)  to  predict  the  varying  ripple  height  q(t)  for  time- 
varying  wave  inputs  is  described  in  Section  6.  DO05  did  not  present  similar  data  for  the 
rate  of  change  of  wavelength,  although  they  remarked  that  the  ripple  length  reached 
equilibrium  at  about  the  same  time  as  the  ripple  height.  On  this  basis,  we  will  take  an 
analogous  equation  for  wavelengths,  X,  to  that  given  for  heights  in  Eq  (4.5),  and  use  the 
same  Eq  (4.6)  for  the  rate  of  change  of  wavelength.  By  analogy,  but  without 
justification  due  to  lack  of  data,  we  will  use  analogous  equations  for  the  variation  of 
ripple  orientation  in  response  to  changes  in  wave  direction. 
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5.  New  current-only  ripple  geometry  predictor 

5.1  NEW  CURRENT-RIPPLE  EQUILIBRIUM  PREDICTOR 

As  shown  in  Section  3,  there  are  far  fewer  current-generated  than  wave-generated  ripple 
predictors,  and  there  have  not  been  definitive  intercomparison  tests  of  them.  We 
therefore  again  developed  a  new  predictor  for  this  project. 

In  developing  a  new  current-ripple  predictor  we  draw  heavily  on  laboratory  experiments 
by  Baas  (1993),  Whitehouse  et  al  (1998)  and  Damgaard  et  al  (2003). 

The  PhD  thesis  of  Baas  (1993)  [B93]  describes  an  extensive  series  of  laboratory 
experiments  in  a  continuously-circulating  race-track  flume  at  the  University  of  Utrecht, 
the  Netherlands.  Two  sediments  were  used:  a  very  fine  sand  with  d50  =  0.095mm,  and  a 
fine  sand  with  d5o  =  0.238mm.  The  experiments  were  run  for  sufficiently  long  durations 
that  both  the  evolution  of  the  ripples  and  their  equilibrium  form  could  be  studied  with 
confidence.  A  wide  range  of  current  speeds  was  used.  A  comprehensive  compilation  of 
data  from  previous  researchers  was  also  made,  and  used  to  assist  interpretation  of  the 
experimental  results. 

Whitehouse  et  al  (1998)  [WMS98]  performed  experiments  in  a  reversing  flume  at 
HR  Wallingford,  UK.  A  single  sediment  was  used  with  dso  -  0.510mm.  Initial 
experiments  were  conducted  with  steady  currents,  with  step  changes  between  current 
speeds,  to  investigate  the  equilibrium  ripple  geometry,  and  the  evolution  from  a  plane 
bed.  Subsequent  experiments  were  performed  with  a  flow  that  reversed  sinusoidally 
with  a  period  of  12.5  hours  for  three  tidal  periods.  This  was  regarded  as  a  direct  full- 
scale  simulation  of  the  bottom  0.1m  or  so  of  a  full-scale  tidal  flow  in  the  sea.  The  aim 
was  to  investigate  tidal  variations  in  ripple  geometry. 

Damgaard  et  al  (2003)  [DSPW03]  performed  experiments  in  a  sloping  sediment  duct  at 
HR  Wallingford,  UK.  Two  sediments  were  used:  a  well-sorted  fine  sand  with  d50  — 
0.237mm  and  d9o/d !0  =  2.4,  and  a  widely-graded  fine  sand  with  d5o  -  0.231mm  and 
d9(>/d  io  =  4.2.  The  main  aim  was  to  investigate  the  effect  of  steep  slopes  (up  to  20°)  on 
sediment  transport  rates,  but  the  ripple  properties  were  also  measured.  Only  the  tests 
with  a  horizontal  bed  were  used  for  the  present  purposes. 

WMS98  and  DSPW03  do  not  appear  in  Tables  1  to  7,  because  they  have  not  been  used 
in  existing  intercomparisons. 

B93  presented  convincing  evidence  for  his  hypothesis  that,  within  the  rippled-bed 
regime  (i.e.  for  all  current  speeds  greater  than  the  threshold  of  motion  but  less  then  the 
start  of  wash-out),  a  given  grain-size  will  develop  ripples  of  a  constant  height  and 
wavelength  provided  that  the  current  speed  is  maintained  for  a  long  enough  time.  This 
is  contrary  to  most  earlier  interpretations,  in  which  ripple  height  was  assumed  to 
increase  with  increasing  current  speed  to  a  certain  point,  then  decrease  due  to  wash-out. 
B93  also  found  that  current  speeds  greater  than  a  certain  value  led  to  a  flat  bed  with 
sheet  flow  (upper  stage  plane  bed)  for  very  fine  sand  with  d50  =  0.095mm,  but  led  to  a 
transition  to  dunes  of  greater  height  and  wavelength  than  the  ripples  for  fine  sand  with 
d50  =  0.238mm.  B93  proposed  expressions,  based  on  his  own  data  and  that  from  other 
sources,  for  the  maximum  ripple  height  and  wavelength. 


r|max  =  3.4  log|0(d50)+  18 


(5.1) 
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^max  ~  75.4  logio('dso)  +  197 


(5.2) 


where  T|max,  7^ax  and  d5o  are  all  in  mm.  These  dimensionally  inhomogeneous  equations 
can  be  re-cast  in  homogeneous  form  if  it  is  assumed  that  (a)  the  length-scale 
determining  r|max  and  is  d50,  (b)  that  the  remaining  influence  of  d50  is  expressed  by 
the  commonly-used  dimensionless  grainsize  D*,  where 


(5.3) 


Equations  (5.1)  and  (5.2)  can  then  be  written  as: 


(5.4) 


(5.5) 


where  In  is  the  natural  logarithm. 

However,  these  expressions  do  not  fit  closely  to  the  measurements  of  WMS98  and 
DSPW03.  Instead,  the  following  expressions  are  proposed  here,  which  fit  closely  to  the 
data  of  B93,  WMS98  and  DSPW03,  and  reasonably  well  to  the  other  data  collated  by 
B93: 

T|max  =  d50 . 202  D.~0554  ,  for  1 .2  <  D.  <  14  (5.6a) 

kmax  =  d50  (500  +  1881  D/15)  ,  for  1.2  <D.  <  14  (5.6b) 

For  D*  >  14  (d50  >  0.7mm)  ripples  are  indistinct  or  do  not  form. 

The  prediction  curves  and  all  the  data  mentioned  above  are  plotted  in  Figures  8  and  9. 
The  equations  were  fitted  exactly  to  the  very-fine-sand  data  of  B93  and  the  coarsc-sand 
data  of  WMS98.  They  fit  quite  closely  to  the  fine-sand  data  of  B93  and  DSPW03,  and 
reasonably  well  to  the  other  data. 

The  effect  of  wash-out  at  large  current  speeds  is  included  by  again  using  the  data  of 
B93,  WMS98  and  DSPW03  together  with  a  plot  by  Van  den  Berg  and  Van  Gelder 
(1989)  reproduced  by  Van  Rijn  (1993).  All  these  are  shown  in  Figure  10,  with  lines 
added  as  an  approximation  to  the  start  of  wash-out  and  the  completion  of  wash-out  and 
the  start  of  sheet-flow  (with  or  without  dunes).  The  lines  are  given  in  terms  of  the  skin- 
friction  Shields  parameter  0'c  as  a  function  of  D*.  For  present  purposes,  0'c  is  defined 
by: 


(5.7a) 


0.40 


(5.7b) 


ln(h/z0)- 1 


z0  =d50/12 


(5.7c) 
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although  the  plot  by  Van  den  Berg  and  Van  Gelder  used  Zo  =  d90/10  (the  differences  are 
not  great).  Here  U  is  depth-averaged  current  speed,  CD  is  drag  coefficient,  and  Zo  is 
grain-related  bed  roughness  length. 

The  proposed  limits  of  wash-out  and  sheet-flow  are: 

e;v0  =  1 .66  D2 3  for  D»  >  1 .58  (5.8a) 

0'f  =  2.26  D(' 3  for  D.  >  1 .58  (5.8b) 

0'wo  =  0.916  and  0'sf  =  1.25  for  D.  <1.58  (5.8c) 


Applying  a  linear  reduction  in  ripple  height  with  0'c,  the  expression  for  equilibrium 
ripple  height  including  wash-out  is: 


r|eq  =  pre-existing  value 

^leq  =  max 


Beq  =  Bmax 
^leq  =0 


e;r  -  e; ' 


for  0  <  0'c<  0cr 
for  0rr  <  0'  <  0\, 


for  0'wo  <  0'  <  0' 


for  0'c  >  0'f 


(5.9a) 

(5.9b) 

(5.9c) 

(5.9d) 


The  wavelength  is  assumed  to  be  unaffected  by  wash-out  (i.e.  only  the  ripple  steepness 
declines)  and  is  given  by: 


^eq  -  ^max 


(5.10) 


No  account  is  taken  here  of  the  lengthening  of  X  on  transition  to  dunes  for  fine  and 
medium  sands.  The  expressions  for  t|max,  Xmax,  0fwo  and  0'r  are  given  by  Eqs  (5.6a), 
(5.6b),  (5.8a,  c)  and  (5.8b,  c)  respectively. 


5.2  TIME-EVOLUTION  OF  CURRENT-RIPPLES 

As  with  wave-ripples,  current-ripples  alter  their  height,  wavelength  and  orientation  with 
changing  current  speed  and  direction.  (They  also  alter  their  asymmetry,  but  that  is  not 
considered  here.)  The  speed  of  ripple  evolution  decreases  as  the  current  speed 
decreases,  and  is  zero  for  currents  below  the  threshold  of  motion  of  sediment  (frozen 
ripples). 

By  analogy  with  the  approach  used  for  wave-ripples  in  Section  4.2  (c.f.  Eq  4.5),  the 
rates  of  change  of  ripple  height  and  wavelength  are  given  by: 


dt  Tc  V  eq  ' 

(5.11) 

dt  Tc  ”  ' 

(5.12) 
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where  r|eq  is  given  by  (5.9a-d),  and  Xm,  is  given  by  Eq  (5.10).  B93  found  that  the  ripple 
height  evolved  faster  than  the  wavelength,  and  hence  two  different  parameters,  Pn  and 
px,  are  employed.  The  basic  time-scale  Tc  is  less  easy  to  identify  for  currents  than  it  was 
for  waves,  where  wave  period  T  was  the  natural  choice.  Here  we  take  Tc  to  be  related  to 
the  time  taken  for  the  volume  of  a  ripple  per  unit  width  (proportional  to  r|ma)t  .  A,rax)  to  be 
delivered  by  the  bedload  transport  rate,  qb.  The  bedload  transport  can  be  written  in  the 
form  (Soulsby,  1997,  p.  158): 


0>  = 


qb 

[g(s-l)d»F 


func  (e;,0cr) 


(5.13) 


The  adaptation  time  in  Eqs  (5.1 1)  and  (5.12)  can  then  be  written  as: 
T  X 

ic  _  '  I  max  '^max 

P* 


[g(s-l)d’0f  'func(e;,0cr) 

from  which  the  time-scale  Tc  and  rate-coefficient  px  (where  x  =  r|  or  X)  can  be  separated 
as: 


qb 

max  max 


X  = 


(5.15) 


and  (3X  -  func(0'Cj  0cr),  where  x  =  r\  or  X,  and  Pn  and  Px  are  different  functions  of  0fC}  0cr. 


B93  tabulated  data  from  his  experiments  for  the  time  taken  for  ripples  to  achieve  99%  of 
their  final  height  or  wavelength.  These  data  have  been  used  here  to  devise  and  calibrate 
expressions  for  pn  and  p*: 


20(e;-eJi 

2.5+(e;-eJ5 

12(8;  -ej5 

2.5+(e;-e„)''s 


(5.16) 


(5.17) 


These  curves  are  shown  in  Figures  1  la  -  d5  together  with  the  data  from  B93,  and  give  a 
reasonably  good  fit. 

In  the  absence  of  data,  the  orientation  (p  of  the  ripples  is  assumed  to  change  at  the  same 
rate  as  the  wavelength.  Thus  the  evolution  is  taken  to  be  given  by: 

f  =  (5.1*) 

where  (peq  is  the  direction  of  the  instantaneous  current. 
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6.  Time-stepping  procedure 

6.1  AIMS  OF  PREDICTOR 

The  aim  is  to  develop  a  ripple  predictor  that  predicts:  time-varying  ripple  heights, 
lengths  and  orientations,  for  any  sandy  sediment,  driven  by  any  time-series  of:  water 
depths,  wave  heights,  periods  and  directions,  and  current  speeds  and  directions.  It  must 
take  account  of:  evolution  (“history”),  threshold-of-motion,  and  wash-out  effects,  as 
well  as  bio-degradation  of  ripples.  It  must  cover  both  wave-generated  and  current- 
generated  ripples,  switching  between  them  (in  an  evolving  sense)  depending  on  which 
forcing  is  dominant. 

6.2  BASIC  EQUATIONS 

The  full  predictor  is  based  on  the  methods  and  equations  derived  in  Sections  4  (for 
waves)  and  5  (for  currents).  It  is  set  out  in  step-by-step  algorithmic  form  in  Appendix 
A. 

The  equilibrium  height  and  wavelength  for  wave-generated  ripples  are  given  by  Eqs 
(4.1)  -  (4.3),  with  the  threshold  constraint  given  after  Eq  (4.3). 

The  equilibrium  height  and  wavelength  for  current-generated  ripples  are  given  by  Eqs 
(5.6)  -  (5.10),  with  the  threshold  constraint  given  by  Eq  (5.9a). 

These  methods  include  both  threshold-of-motion  and  wash-out  effects.  The  differential 
equations  expressing  the  evolution  of  ripple  height  r|(t),  length  A,(t)  and  orientation  cp(t) 
are  given  by  Eq  (4.5)  for  wave-generated  ripple  height,  with  analogous  equations  for 
length  and  orientation;  and  by  Eqs  (5.11),  (5.12)  and  (5.18)  respectively  for  current- 
generated  ripples. 

The  rate-of-change  parameters  in  these  equations  are:  p/T  for  wave-generated  ripple 
height,  length  and  orientation;  and  Pt/Tc,  p^/Tc  and  px/Tc  respectively  for  current- 
generated  ripple  height,  length  and  orientation.  In  these: 

•  p  is  given  by  Eq  (4.6),  where  \| /  is  given  by  Eq  (3.3) 

•  T  is  wave  period  (see  discussion  in  Section  3.4  on  selection  of  T) 

•  pn  is  given  by  Eq  (5.16),  where  0'c  is  given  by  Eqs  (5.7a)  to  (5.7c) 

•  p>.  is  given  by  Eq  (5.17) 

•  Tc  is  given  by  Eq  (5.15),  where  r|max  and  X^  are  given  by  Eqs  (5.6a)  and  (5.6b). 

The  differential  equations  for  ripple  orientation  cp  (=  direction  of  normal  to  ripple  crests) 
are  an  interim  approximation  only.  They  do  not  behave  correctly  when  (p  passes 
anticlockwise  through  0°N;  for  example,  if  cpeq  makes  a  step  change  from  10°N  to 
350°N,  cp(t)  rotates  clockwise  via  South  instead  of  anticlockwise  through  North.  An 
alternative  approach,  expressing  wavelength  X  as  a  vector  ( X Xy)  was  not  an 
improvement.  In  fact,  any  such  approach  with  a  single  ripple  train  is  unrealistic, 
because  a  spatially  extended  ripple  field  cannot  swing  around  as  a  whole  when  the  wave 
(or  current)  direction  changes.  The  present  approach  should  therefore  be  used  with 
caution,  and  preferably  only  for  wave  (or  current)  directions  within  a  180°  sector  that 
excludes  0°.  It  is  hoped  to  devise  a  more  robust  method  in  the  future. 
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6.3  WAVES  VERSUS  CURRENT  RIPPLES 

The  present  version  of  the  time-stepping  procedure  adopts  a  bi-polar  approach  to  the 
relative  effects  of  waves  and  currents  in  generating  ripples.  This  amounts  to  a  switch 
from  the  wave-only  to  the  current-only  expressions,  and  vice  versa,  depending  on  the 
dominant  value  of  Shields  parameter.  A  similar  approach  was  used  by  Lyne  et  al  (1990) 
based  on  the  ratio  of  current  speed  to  wave  orbital  velocity.  However,  the  present 
approach  is  proposed  on  the  stronger  basis  that  bedload  sediment  transport  rates  (which 
are  the  driving  mechanism  for  ripple  orientation)  depend  directly  on  Shields  parameter. 

Thus,  comparing  0'w  for  waves  (Section  4)  with  0'c  for  currents  (Section  5): 

•  if  e\v  >  0'c,  the  wave-generated  expressions  summarised  in  Section  6.2  are  used 

•  if  0'c  >  0'w,  the  current-generated  expressions  summarised  in  Section  6.2  are  used 
where  0'wis  given  by  Eq  (3.4)  and  0'c  by  Eq  (5.7a,  b,  c). 

Li  and  Amos  (1998)  proposed  expressions  for  ripple  geometry  due  to  combined  waves 
and  currents,  and  gave  a  three-fold  categorisation  of  the  relative  effects  based  on  the 
ratio  0'w/0'c.  Neither  Lyne  et  al  (1990)  nor  Li  and  Amos  (1998)  used  a  time-evolving 
predictor;  they  just  assumed  that  the  ripple  geometry  was  always  in  equilibrium  with  the 
instantaneous  hydrodynamic  conditions. 

6.4  BIO-DEGRADATION  OF  RIPPLES 

Marine  creatures  burrowing  in  the  sea-bed  sediments  can  cause  a  ripple  field  to  flatten. 
This  is  especially  noticeable  at  slack  water  for  (tidal)  current-generated  ripples,  or  under 
calm  conditions  for  wave-generated  ripples.  Amos  ct  al  (1988)  observed  at  a  depth  of 
22m  on  the  Scotian  shelf  that  re-working  by  burrowing  organisms  could  mix  the  top 
10cm  or  so  of  the  sea-bed  within  4  to  6  hours,  thereby  flattening  current- formed  ripple 
patterns  in  periods  when  both  tidal  flow  and  wave  activity  were  sub-threshold.  Baas 
(1993)  described  video  observations  of  the  bottom  of  a  subtidal  channel  in  the  Dutch 
Wadden  Sea  in  which  intensive  bioturbation  obliterated  ripples  formed  during  one  tide. 
He  also  reported  geological  evidence  of  current  ripples  in  turbidite  beds  found  in  the 
Pyrenees  Mountains  (Spain),  many  of  which  showed  evidence  of  bioturbation.  In  the 
SAX04  experiment,  fish  were  observed  to  pock-mark  the  sea-bed  sediments,  disturbing 
the  ripple  pattern.  Wheatcroft  (1994)  observed  the  bottom  roughness  of  a  silt  bed  in 
90m  of  water  on  the  central  California  shelf.  The  rms  heights  of  the  bottom  roughness 
were  found  to  vary  over  short  periods  of  time  (12  hours).  In  general,  all  physical 
bedforms  could  be  destroyed  by  bioturbation  processes  in  periods  of  hours  to  days. 
Marine  creatures  can  also  create  bed  topography.  For  example,  the  dominant  bed 
roughness  in  the  intertidal  Eden  estuary  in  eastern  Scotland  is  caused  by  a  pattern  of  the 
eroded  worm-casts  of  burrowing  marine  worms  (personal  observation). 

The  detailed  effects  of  bioturbation  depend  on  the  organism  involved.  However,  a 
simple  method  is  proposed  here  for  simulating  the  effects  of  both  ripple  flattening  and 
residual  biologically-induced  bed  roughness.  The  reduction  in  height  of  ripples  is 
expressed  by  a  biological  half-life  T /hb  representing  the  time  taken  for  biological  activity 
to  reduce  the  ripple  height  to  one  half  of  its  initial  value  in  conditions  when  waves  and 
currents  are  below  the  threshold  of  motion.  The  residual  bio-roughness  is  expressed  by 
T]b,  representing  the  average  “trough  to  crest”  height  of  the  biologically  induced  bed 
features. 

These  are  implemented  as  follows: 
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1. 

2. 


3. 

4. 


A  switch  Swb  is  given  the  value  1  if  bio-degradation  is  included,  and  0  if  it  is  not 
Equation  (4.5)  is  modified  to: 


*!  =  £(,,  n)_JL.Sw 

dt  TVLq  J  T. 


(6.1) 


where  Tb  =  T./3>b/ln(2)  =  1.443  T»/j>b 

The  ripple  height  r|(t)  is  set  to  the  maximum  of  qb  or  the  value  derived  from 
integration  of  Eq  (6.1). 

The  same  rules  are  applied  to  Eq  (5.1 1)  if  the  ripples  are  current-generated. 


The  biological  term  in  Eq  (6.1)  causes  the  ripple  height  to  decrease  exponentially  with 
time.  The  ripple  length  and  orientation  are  not  modified  by  biological  effects  in  the 
present  treatment. 


The  choice  of  values  for  T>/Jjb  and  r|b  must  be  made  using  knowledge  about  the  biological 
activity  of  the  study  site.  At  present,  it  is  not  easy  to  recommend  values,  but  ultimately 
this  might  be  possible  by  examining  evidence  from  a  number  of  representative  sites. 


6.5  SOLUTION  METHOD 

The  various  differential  equations  for  the  evolution  of  ripple  height,  wavelength  and 
orientation,  produced  by  waves  or  by  currents,  can  all  be  written  in  the  general  form: 

A7  =  a(t)—  b(t)-  x(t)  (6.2) 

dt 

where  x  represents  either  height,  wavelength  or  orientation.  The  time-varying 
coefficients  a(t)  and  b(t)  are  related  to  the  time-varying  inputs  of  wave  heights  and 
periods,  water  depths,  and  current  speeds.  For  example,  taking  Eq  (6.1)  where  x  =  q,  the 
coefficients  arc: 


a 


6)  = 


(6.3) 


b(0 


+  —  ■  Sw 
T 


(6.4) 


and  p,  T  and  qeq  are  all  time-varying  with  the  wave  conditions.  The  observed  values  of 
the  wave  and  current  conditions  are  usually  recorded  at  discrete  intervals  of  time.  At,  in 
lab  and  field  measurements,  where  At  may  or  may  not  be  constant  throughout  the 
experiment.  The  simplest  approach  to  implementing  Eq  (6.2)  is  to  treat  the  coefficients  a 
and  b  as  being  fixed  at  their  value  a(t)  and  b(t)  over  the  time  interval  t  to  t  +  At.  Eq  (6.2) 
then  has  the  analytical  solution: 


^  ill 


X;  + 


vb> 


(l-exp(-b,At)) 


(6.5) 


where  x;  =  x(t),  Xh  =  x(t  +  At),  a;  =  a(t)  and  b\  -  b(t).  Similar  approaches  to  Eq  (6.5) 
were  used  by  Baas  (1993)  for  current-induced  ripples,  Wilbers  (2004)  for  sub-aqueous 
dunes,  and  ABP  (2004)  for  ripples  under  random  waves.  However,  this  approach 
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neglects  the  changes  in  a  and  b  during  the  time  interval  At,  and  hence  is  not  accurate  if  a 
or  b  changes  significantly  in  time  At. 

The  improved  approach  adopted  here  is  to  solve  Eq  (6.2)  for  a  and  b  varying  linearly 
with  time  from  as  to  ai+i,  and  bi  to  bj+i.  This  does  not,  in  the  general  case,  have  an 
analytical  integral  analogous  to  Eq  (6.5).  Various  solution  schemes  were  tried, 
including  an  explicit  forward-stepping  (Euler)  scheme,  two  implicit  schemes,  and  a 
standard  Runge-Kutta  integration  scheme.  AH  of  these  became  unstable  if  (bAt)  became 
too  large.  The  solution  adopted  was  a  Runge-Kutta  integration  with  an  “over-ride”  if 
bAt  >  1.6.  In  the  latter  case,  the  analytical  solution  Eq  (6.5)  is  applied.  However,  it  was 
found  in  the  tests  with  lab  and  field  data-sets  described  in  Section  7  that  the  “over-ride” 
only  takes  effect  for  a  very  small  proportion  of  the  time.  The  adopted  scheme  (Runge- 
Kutta  with  over-ride)  appears  to  be  stable  for  all  inputs  tested.  Details  are  given  in 
algorithmic  form  in  Appendix  A. 

The  time-stepping  procedure  thus  has  the  following  elements  at  each  point  in  time  t,: 

•  calculate  equilibrium  ripple  height,  wavelength  and  orientation  due  to  wave 
forcing,  together  with  P  and  0'w 

•  calculate  equilibrium  ripple  height,  wavelength  and  orientation  due  to  current 
forcing,  together  with  Pn,  P*.  and  0'c 

•  decide  wave  or  current  dominance,  depending  on  0'w  >  or  <  0'c 

•  set  values  of  coefficients  of  aj  and  bj  for  time-step  f,  based  on  wave  or  current 
dominance 

•  calculate  Xi+i  from  Eq  (6.2)  using  Runge-Kutta  integration  (with  over-ride  using  Eq 
(6.5)  if  bAt  >  1 .6)  for  each  of  ripple  height,  wavelength  and  orientation 

•  move  to  next  time-step  (tj+i). 

To  test  the  algorithm,  and  act  as  a  demonstrator,  the  above  procedure  (as  detailed  in 
Appendix  A)  has  been  implemented  in  an  Excel  spreadsheet.  It  would  be  a 
straightforward  matter  to  implement  the  procedure  as  a  FORTRAN  subroutine  (or  other 
preferred  computer  language).  The  Excel  spreadsheet  has  the  title  “Ripple  evolution 
waves  and  currents  V2.0”.  The  first  worksheet  gives  instructions  on  its  use.  Examples 
of  the  input  sheets  are  given  in  Appendix  B.  Because  field  and  lab  data-sets  are  so  very 
varied  (different  recording  intervals,  durations,  measures  of  wave  heights,  period  or 
orbital  velocity,  depth-averaged  or  point  current  measurements)  it  is  expected  that  a  user 
will  need  to  tailor  the  format  of  the  inputs,  and  the  scales  etc.  of  output  plots,  to  specific 
applications.  It  is  research-level,  not  operational,  software. 
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7.  Tests  of  Mark  2  ripple  predictor  against  data 

7.1  TEST  WITH  WAVE-RIPPLE  EVOLUTION  DATA 

The  time-evolving  wave-ripple  predictor  has  been  tested  against  the  detailed  time-series 
of  ripple  height  and  wavelength  in  response  to  step-changes  in  wave  conditions 
measured  by  DO05  (kindly  provided  by  Drs  Doucette  and  O’Donoghue).  In  the 
predictor,  the  values  of  “wave  height”  used  as  input  were  selected  to  match  the  orbital 
diameters  (D  =  2A)  produced  in  the  Oscillating  Water  Tunnel  and  tabulated  by  DO05. 
The  periods  were  set  to  the  values  tabulated  by  DO05,  and  the  sediment  d50  was  as 
reported  by  DO05.  All  the  experiments  used  in  these  tests  were  for  regular  asymmetric 
waves. 

Figure  12  shows  an  experiment  with  ripple  growth,  followed  by  decay  under  reduced 
wave  conditions.  The  bed  was  initially  flat,  and  oscillatory  velocities  with  period  T  =  5s 
and  D  =  0.99m  were  imposed  for  40  minutes.  At  this  time  the  velocities  were  changed 
to  T  =  4s  and  D  =  0.48m,  and  imposed  for  a  further  77  minutes.  The  predictor  simulates 
the  rate  of  growth  in  ripple  height  and  wavelength  reasonably  well,  but  slightly 
underestimates  the  ultimate  (equilibrium)  values.  It  also  simulates  the  rate  of  decrease 
reasonably  well  after  the  change  in  inputs,  this  time  slightly  overestimating  the  final 
values. 

Figure  13  shows  an  experiment  with  ripple  growth  followed  by  further  growth  under 
increased  wave  conditions.  The  bed  was  initially  flat,  and  oscillatory  velocities  with 
T  -  2s  and  D  =  0.24m  were  imposed  for  42  minutes,  after  which  time  the  ripples  had 
attained  equilibrium.  At  this  time  the  velocities  were  changed  to  T  =  4s  and  D  =  0.48m, 
and  imposed  for  a  further  67  minutes.  The  rather  weak  wave  action  at  the  start  did  not 
start  to  form  vortex  ripples  in  the  flattened  bed  for  the  first  15  minutes.  The  ripple 
predictor  was  “switched  off’  during  this  period  when  rolling-grain  ripples  were  present 
(a  common  feature  at  the  start  of  laboratory  experiments  with  flattened  beds).  The 
predictor  simulates  the  rate  of  growth  reasonably  well,  but  slightly  over-prcdicts  the 
ultimate  height  and  wavelength  after  42  minutes.  The  predictor  also  simulates  the 
renewed  growth  after  the  change  in  wave  conditions  reasonably  well,  and  predicts  both 
the  final  height  and  wavelength  reasonably  accurately. 

Figure  14  shows  three  separate  experiments  separated  on  the  plots  by  gaps  (for  ease  of 
plotting,  these  are  not  consecutive  experiments).  The  first  shows  growth  from  an 
initially  flat  bed  for  wave  inputs  of  T  =  6s  and  D  -  0.57m.  The  growth  is  very  slow 
under  these  weak  wave  conditions,  starting  with  about  40  minutes  before  vortex  ripples 
begin  to  form  on  the  flattened  bed,  during  which  time  the  ripple  predictor  was  “switched 
off’.  After  40  minutes,  the  predictor  accurately  simulates  both  the  growth  rates  and  the 
final  height  and  wavelength  over  the  next  160  minutes.  The  second  experiment  is  the 
one  used  in  the  latter  part  of  Figure  13,  but  with  the  starting  ripple  height  and 
wavelength  set  in  the  predictor  to  the  observed  initial  values.  Likewise,  the  third 
experiment  is  the  same  as  the  latter  part  of  Figure  12,  with  the  predictor  started  from 
observed  values.  The  correspondence  of  the  predictions  and  observations  are  improved 
in  these  cases  compared  with  Figures  13  and  12,  where  the  predictor  used  its  own 
predicted  initial  values  from  the  end  of  the  prior  experiments. 

The  performance  of  the  time-evolving  wave-ripple  predictor  is  thus  reasonably  good.  It 
predicts  the  rate  of  change  quite  accurately  (although  the  tests  are  for  the  data-set 
against  which  it  was  calibrated).  It  also  gives  acceptably  good  predictions  of  the 
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equilibrium  ripple  height  and  wavelength.  Some  of  the  deviations  may  be  due  to  the 
orbital  velocities  in  the  experiments  being  strongly  asymmetrical,  which  the  predictor 
was  not  specifically  designed  for. 

7.2  TEST  WITH  CURRENT-RIPPLE  “TIDAL”  DATA 

The  expressions  for  ripple  evolution  were  derived  in  Section  5.2  using  only  the  data  of 
Baas  (1993)  for  the  growth  of  ripples  from  an  initially  flat  bed.  It  has  been  tested 
against  the  more  realistic  case  of  the  evolution  of  ripples  under  continuously  changing 
currents,  using  an  independent  data-set  from  the  tidally  reversing  flume  experiments  of 
Whitehouse  et  al  (1998).  These  started  with  an  initially  flat  bed,  and  measured  the 
ripple  development  over  three  tidal  cycles  (37.5  hours).  The  six  half-cycles  are  denoted 
as  1  to  6  in  the  commentary  on  the  results  given  below.  Here  we  have  used  the  current 
speed  measured  at  a  height  of  0.10m  in  the  flume  to  drive  the  current-induced  ripple 
predictor.  The  equilibrium  height  and  wavelength  are  given  by  Eqs  (5.9a-d)  and  (5.10). 
The  time  evolution  is  obtained  by  a  time-stepping  integration  of  Eqs  (5.1 1)  and  (5.12), 
in  which  Tc  is  given  by  Eq  (5.15),  and  Pn  and  (3^by  Eqs  (5.16)  and  (5.17)  respectively. 

It  is  seen  in  Figure  15  that  the  time-evolving  ripple  predictor  manages  to  predict  many 
of  the  features  seen  in  the  data:  the  growth  from  a  flat  bed,  the  effects  of  a  slight 
asymmetry  in  the  tidal  cycle,  and  the  partial  wash-out  that  occurs  at  the  peak  of  the 
negative  (stronger)  half-cycles.  The  following  features  can  be  seen.  In  the  first  half¬ 
cycle,  the  model  starts  to  form  ripples  from  an  initially  flat  bed  at  about  the  right  time. 
The  predicted  growth  rate  for  ripple  height  is  a  bit  too  fast,  but  is  about  right  for  ripple 
length.  The  slight  washout  at  maximum  velocity,  which  briefly  halts  the  growth,  is 
predicted  correctly.  The  slack-water  dwell  between  the  first  and  second  half-cycles 
(when  velocity  is  below  threshold)  occurs  for  about  the  right  time,  although  the  height 
and  length  are  slightly  over-estimated.  The  similar  dwells  at  subsequent  slack  waters 
are  predicted  quite  accurately  for  both  height  and  length.  After  the  first  half-cycle,  the 
predicted  ripple  height  continues  to  grow  slowly,  in  line  with  the  observations,  attaining 
almost  the  equilibrium  value  towards  the  end  of  the  second,  fourth  and  sixth  half-cycles. 
The  ripples  start  to  wash  out  near  the  maximum  velocity  in  half-cycles  1,  3  and  5, 
because  the  velocity  is  slightly  larger  than  in  half-cycles  2,  4  and  6.  The  model 
correctly  predicts  the  reductions  in  height,  but  does  not  predict  the  increases  in 
wavelength  at  these  times  (this  is  probably  the  start  of  a  transition  to  dunes,  which  is  not 
included  in  the  model).  At  25  hours  the  flume  was  switched  off  for  over  a  month, 
during  which  time  the  ripple  height  slumped  by  about  6mm  due  to  artificial  disturbance. 
The  standard  run  of  the  model  does  not  (of  course)  reflect  this,  but  in  a  test  in  which  the 
height  was  re-set  at  25.5  hours  to  the  slumped  value,  the  model  re-built  the  height  very 
much  as  observed.  Thus  the  model  reproduces  most  of  the  features  seen  in  the 
observations.  The  predicted  ripple  heights  are  nearly  always  within  5mm  of  the 
observed  values,  and  the  wavelength  generally  within  about  5cm.  The  predictor  does 
not  capture  the  increase  in  wavelength  at  peak  flows,  as  this  is  associated  with  dune 
development,  which  the  predictor  is  not  designed  to  include. 

Overall,  the  performance  of  the  predictor  is  seen  to  be  good.  In  fact,  it  even  elucidates 
aspects  of  the  measurements  that  were  not  apparent  from  the  data  alone. 

7.3  TEST  WITH  WAVE-PLUS-CURRENT  SYNTHETIC  DATA 

A  test  was  made  of  the  full  wave -plus-current  predictor  to  check  that  it  can  predict  the 
features  intended,  such  as  ripple  growth  and  decay,  wash-out,  frozen  ripples,  and 
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switching  from  wave-dominated  to  current-dominated  conditions  (and  vice  versa).  A 
time-series  of  waves  and  currents  was  synthesised  in  which: 

•  Sediment  d50  =  0.2mm 

•  Water  depth  =  10m  (constant) 

•  Wave  heights  vary  sinusoidally  between  0m  (calm)  and  2m  (storm)  with  a  period 
of  120  hours  (5  days)  between  storms 

•  Wave  period  varies  sinusoidally  between  4s  and  8s  in  line  with  heights 

•  Wave  direction  varies  sinusoidally  between  80°N  and  140°N  in  line  with  heights 

•  Current  speed  varies  sinusoidally  with  amplitude  0.8ms'1  and  reverses  rcctilinearly 
between  90°N  and  270°N,  with  period  12  hours  (tidal) 

•  Bio-degradation  switched  OFF 

•  Initial  values  of  ripple  height,  wavelength  and  orientation  are  set  to  0.02m,  0.2m 
and  100°N. 

The  resulting  ripple  evolution  repeats  with  a  5-day  period  (Fig  16).  The  following 
features  can  be  seen.  After  an  initial  rapid  transient  at  the  start  of  Day  0,  ripples  of 
height  0.017m,  wavelength  0.13m  and  tidally  oscillating  orientation  are  formed  under 
the  current-dominant  conditions  while  the  waves  are  weak.  At  the  peak  velocity  of  each 
tidal  half-cycle,  the  ripple  height  reduces  to  0.008m  due  to  partial  wash-out,  then 
recovers  back  to  0.017m  as  the  current  speed  drops.  This  behaviour  repeats  for  four 
tidal  half-cycles  until  the  wave  heights  pick  up  at  the  start  of  Day  1. 

From  Day  1.0  to  Day  1.5,  there  are  alternating  periods  of  wave  dominance  and  current 
dominance,  with  the  ripple  height,  wavelength  and  orientation  all  varying  accordingly. 
At  about  Day  1.5  the  waves  assert  their  dominance,  creating  ripples  of  height  0.032m, 
wavelength  0.22m  and  orientation  110°N.  As  the  waves  become  larger  (increasing 
storm),  the  ripples  progressively  start  to  wash  out,  reducing  the  height  to  0.015m  and  the 
wavelength  to  0.15m.  At  the  same  time,  the  orientation  of  the  ripples  swings  round 
from  110°N  to  140°N,  following  the  wave  direction.  As  the  waves  then  decrease 
(waning  storm),  the  ripple  height  and  wavelength  recover  to  0.032m  and  0.22m,  and  the 
orientation  swings  back  from  140°N  to  110°N,  again  following  the  wave  direction.  At 
Day  3.5,  another  period  of  alternating  wave  and  current  dominance  ensues,  until  from 
Day  4.0  to  5.0  the  current  takes  charge  and  the  ripples  behave  in  the  same  way  as  at  the 
start.  From  Day  5  to  Day  10  the  pattern  repeats. 

All  the  changes  take  place  relatively  quickly,  because  the  maximum  wave  and  current 
conditions  are  both  strong.  Tests  with  other  values  of  the  inputs  show  different 
behaviour  (not  illustrated).  For  example,  reducing  the  peak  wave  height  to  0.5m  makes 
the  ripples  mainly  current-dominated,  with  waves  only  taking  over  around  their 
maximum  heights.  Reducing  the  peak  current  to  0.5ms'1  (with  maximum  H  =  2m)  gives 
a  largely  wave-dominated  picture  with  no  current-induced  wash-out.  Setting  peak 
current  to  0.41ms‘l  and  peak  wave  height  to  0.5m  (both  just  above  threshold)  gives  a 
very  slowly  changing  pattern,  sometimes  current-dominated  and  sometimes  wave- 
dominated,  and  not  fully  established  by  the  end  of  5  days,  so  the  second  5-day  cycle 
looks  very  different  from  the  first.  All  the  tests  resulted  in  the  behaviour  expected. 

7.4  TEST  WITH  WAVE-PLUS-CURRENT  FIELD  DATA 

It  was  originally  intended  that  the  full  Mark  2  Ripple  Predictor  would  be  tested,  and  if 
necessary  re-calibrated,  against  data  obtained  from  the  ONR  field  programs  SAX99  and 
SAX04.  It  subsequently  transpired  that  the  SAX99  data  comprised  only  “frozen 
ripples”,  for  which  the  wave  and  current  conditions  at  the  time  of  their  formation 
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occurred  before  the  measurement  period.  Although  this  demonstrates  the  importance  of 
ripple  predictors  including  the  “frozen  ripples”  effect,  it  was  not  possible  to  test  the 
predictor  without  knowing  the  foregoing  conditions  The  SAX04  data  appear  to  be 
eminently  suitable  for  test  purposes,  but  were  not  at  a  stage  to  be  made  available  within 
the  present  contract. 

Instead,  as  a  substitute,  a  set  of  measurements  at  hourly  intervals  for  1 1  months  of 
waves,  currents  and  water  depths  recorded  on  a  sea-bed  instrumented  frame  at  a  site  off 
Teignmouth  on  the  south  coast  of  England  has  been  used.  This  was  measured  by 
HR  Wallingford  as  part  of  the  EU  collaborative  research  project  COAST3D,  and  is 
described  fully  in  Whitehouse  (2004).  However,  the  ripple  geometries  at  the  recording 
station  were  not  measured,  so  the  data  can  only  be  used  to  test  whether  the  predictor 
gives  plausible  results.  An  approximately  two-month  section  of  data  from  8  July  to  9 
September  1999  has  been  selected  as  a  test-case,  which  includes  two  storms  with  calm 
periods  in  between,  and  four  spring-neap  cycles  of  water  levels  and  tidal  currents.  The 
inputs  are: 

•  Sediment  d50  =  0. 145mm 

•  Water  depth  varying  tidally  (and  with  spring-neap  variation)  between  about  3m 
and  7m 

•  Wave  height  (H|/io)  reaching  about  2m  in  two  storms  each  of  duration  about  5 
days,  and  reaching  0.5m  to  1.1m  in  four  minor  storms,  with  Hj/io  <  0.5m  in  the  rest 
of  the  time 

•  Wave  period  (T i/io)  generally  between  4s  and  10s,  but  occasionally  reaching  14s 

•  Wave  direction  from  offshore  sector,  generally  between  East  and  South-east  (90°N 
to  135°N) 

•  Current  speed  varying  tidally  up  to  about  0.3m.s_l,  occasionally  reaching  0.4m. s’1 
when  wind-assisted 

•  Current  direction  generally  rotating  clockwise 

•  There  is  a  break  at  Day  227  for  servicing  of  instruments 

•  Days  are  counted  from  l  January  1999. 

Bio-degradation  was  switched  OFF  in  the  predictions,  and  initial  ripple  heights,  length 
and  orientation  were  set  to  0.02m,  0.2m  and  100°N.  The  resulting  predictions  (Figure 
17)  demonstrate  a  number  of  features.  The  ripples  are  almost  entirely  wave-generated. 
On  the  few  occasions  that  the  current  dominates,  it  is  so  weak  that  it  has  little  or  no 
effect  on  the  pre-existing  ripples.  Small  waves  in  days  189  to  206  produce  ripples  with 
heights  between  0.014  and  0.023m,  wavelengths  between  0.09  and  0.16m  and 
orientations  between  90°  and  150°N.  There  are  long  periods  of  frozen  ripples,  when 
neither  the  waves  nor  the  current  exceed  the  threshold  of  motion.  The  first  storm, 
starting  day  206,  initially  builds  ripples,  but  then  partially  washes  them  out,  with  some 
periods  of  recovery  near  low  water  when  the  wave  heights  are  reduced.  As  the  storm 
wanes,  the  ripple  height  and  length  grow  to  about  0.023m  and  0.15m  respectively. 
From  days  210  to  232  the  ripples  respond  to  minor  storms.  The  second  major  storm 
(days  233  to  237)  partially  washes  out  the  ripples,  which  recover  by  day  235,  and  the 
ripples  re-grow  during  the  waning  storm.  From  days  239  to  245  neither  the  waves  nor 
the  currents  are  strong  enough  to  alter  the  ripples,  which  remain  frozen.  Minor  storms 
in  the  last  6  days  re-activate  the  ripples.  The  main  features  expected  of  wave-generated 
ripples  thus  appear  to  be  (at  least  qualitatively)  predicted  correctly. 

Further  tests  have  been  made  by  scaling  the  measured  Teignmouth  inputs. 
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Stronger  currents  and  weaker  waves  (Figure  18).  The  depths  were  doubled  (weakening 
the  wave  orbital  velocities)  and  the  current  speeds  were  also  doubled.  Most  of  the  time 
the  ripples  are  current-generated  (and  without  any  wash-out),  giving  constant  height  = 
0.015m,  constant  length  =  0.12m,  and  alternating  orientations.  During  the  two  major 
storms,  and  to  a  lesser  extent  the  minor  storms,  the  ripple  height,  length  and  orientation 
become  wave-related. 

Biodegradation  switched  ON  (Figure  19).  The  biodegradation  half-life  Thy,  was  set  to 
50  hours  (corresponding  to  Tb  =  72  hours),  and  the  residual  bio-roughness  was  set  to 
5mm.  Whenever  the  physical  ripple  growth  is  small  or  zero,  the  pre-existing  ripple 
height  decays  due  to  biological  activity.  The  “ripple  height”  does  not  decrease  below 
5mm,  because  of  the  residual  bio-roughness.  The  predicted  length  and  orientation  are 
not  affected  in  these  periods. 
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8.  Conclusions 

A  time-evolving  sea-bed  ripple  predictor  has  been  developed,  based  on  a  mixture  of 
European  and  North  American  data  and  understanding.  It  predicts  the  time-vaiying 
ripple  heights,  lengths  and  orientations  in  sandy  sediments  under  the  action  of  waves 
and  currents.  It  takes  account  of  evolution  (history  effect),  threshold-of-motion  effects 
(frozen  ripples),  wash-out  effects  (leading  to  flat  bed),  and  biological  effects  (bio¬ 
degradation). 

Detailed  achievements  are: 

•  Digest  of  source  publications  (Report  TR  150,  Soulsby  and  Whitehouse,  2005) 

•  Data-base  of  ripple  geometries  and  driving  conditions  (Report  TR  150) 

•  Statistical  analyses  of  ripple  heights  and  wavelengths  (Report  TR  150) 

•  Made  an  “Intercomparison  of  intercomparisons”  -  24  intercomparison  papers  were 
analysed,  which  overall  compared  28  ripple  prediction  methods  (for  wave  and/or 
current-generated  ripples)  with  (subsets  of)  84  data-sets.  There  is  a  need  to  make  a 
definitive  test  of  all  predictors  against  all  data 

•  Devised  Mark  1  Ripple  Predictor  for  equilibrium  ripples  under  (a)  waves,  (b) 
currents 

•  Devised  Mark  2  Ripple  Predictor  for  the  time-evolution  of  ripples 

•  Tested  Mark  2  Ripple  Predictor  against  laboratory'  and  field  data-sets. 

The  Mark  2  Ripple  Predictor  gave  good  agreement  with  quasi-full-scale  laboratory 
measurements  of  ripple  growth  and  decay  due  to  (a)  changing  wave  conditions,  (b) 
changing  current  conditions.  The  tests  against  field  data  were  only  qualitative,  and  tests 
against  field  data  with  time-series  of  ripple  measurements  were  not  possible.  Tests 
against  SAX04  data  await  the  release  of  that  data  in  the  future.  The  data  used  to  derive 
and  test  the  predictor  covered  sediments  in  the  range  of  sizes  0.095mm  <  d5o  <  0.51mm. 
The  predictor  is,  by  extension,  expected  to  give  realistic  results  in  the  grain-size  range 
0.062mm  to  0.7mm.  Finer  sediments  will  increasingly  feel  the  effects  of  electro¬ 
chemical  cohesion,  which  will  limit  ripple  formation.  Sediments  coarser  than  about 
0.7mm  do  not  usually  form  ripples.  No  limits  need  to  be  imposed  on  water  depths  (on 
the  continental  shelf),  wave  conditions  or  current  conditions,  since  the  predictor  is 
designed  to  cope  with  all  values.  The  method  used  for  predicting  ripple  orientation  is  an 
interim  approach,  limited  in  general  to  wave  and  current  directions  that  do  not  pass 
anticlockwise  through  0°N. 

The  Mark  2  Ripple  Predictor  is  presented  as  a  step-by-step  algorithm  in  Appendix  A, 
suitable  for  coding  as  a  computer  subroutine  or  program  (not  done  in  present  work)  in 
any  preferred  language  and  for  interfacing  to  other  (e.g.  acoustic  propagation)  software 
packages. 

An  Excel  spreadsheet  has  been  created  to  test  and  demonstrate  the  predictor.  It  can  be 
tailored  to  data-sets  in  different  forms  to  the  ones  tested,  but  is  not  intended  as  an 
operational  tool.  Individual  results  can  be  fed  into  the  Excel  ripple  visualiser  developed 
in  Phase  1  (Report  TR  150). 
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Table  1  Ripple  prediction  methods  and  intercomparison  papers 


Predictor/Intercomparison 

Code 

Wave  (W), 
Current  (C)  or 
both  (W+C) 

Equilibrium 
(E)  or 

transient  (T) 

Predictor 

Inter¬ 

comparison 

Miller  &  Komar  1980 

MK80 

W 

E 

P 

Nielsen  1981 

N81 

W 

E 

P 

Grant  Sc  Madsen  1982 

GM82 

w 

E 

P 

Clifton  Sc  Dingier  1984 

CD  84 

w 

E 

P 

Tanaka  Sc  Shuto  1984 

TS84 

W+C 

E 

P 

I 

Vongvisessomjai  1984 

V84 

w 

E 

P 

I 

Yalin  1985 

Y85 

c 

E 

P 

I 

McLean  &  Smith  1986 

MS86 

c 

E 

P 

Vongvisessomjai  1987 

V87 

w 

T 

P 

Nelson  &  Smith  1989 

NS89 

c 

E 

p 

Van  Rijn  1989 

VR89 

w 

E 

P 

Lyne  et  al  1990 

L90 

W+C 

E 

P 

1 

Nielsen  1992 

N92 

w 

E 

P 

1 

Wiberg  &  Nelson  1992 

WN92 

c 

E 

P 

I 

Baas  1993 

B93 

c 

E 

P 

I 

Madsen  1993 

Ma93 

w 

E 

P 

Mogridge  et  al  1993 

M93 

w 

E 

P 

1 

van  Rijn  1993 

VRN93 

w 

E 

1 

van  Rijn  et  al  1993 

VR93 

w 

E 

P 

Mogridge  et  al  1994 

M94 

w 

E 

P 

I 

Wiberg  Sc  Harris  1994 

WH94 

w 

E 

P 

I 

Li  et  al  1996 

L96 

W+C 

E 

P 

Tanaka  Sc  Dang  1996 

TD96 

W+C 

E 

P 

Li  &  Amos  1998 

LA98 

W+C 

E 

I 

Traykovski  et  al  1999 

T99 

w 

E 

J 

Khelifa  &  Ouellet  2000 

KO00 

W+C 

E 

P 

1 

Hanes  et  al  200 1 

H01 

W 

E 

1 

O'Donoghue  Sc  Clubb  2001 

ODC01 

W 

E 

I 

Faraci  Sc  Foti  2002 

FF02 

w 

E 

P 

1 

Grasmeijer  2002 

G02 

w 

E 

P 

Foti  &  Faraci  2003 

FF03 

w 

E 

I 

ABP  2004 

ABP04 

w 

T 

P 

1 

Grasmeijer  &  Kleinhans  2004 

GK04 

w 

E 

1 

Williams  et  al  2004 

W04 

w 

E 

I 

Doucette  Sc  O'Donoghue  2005a, b 

DO05 

w 

T 

P 

I 

O'Donoghue  et  al  2005 

OD05 

w 

E 

I 

Smith  &  Sleath  2005 

SS05 

w 

T 

P 

van  Rijn  2005 

vR05 

w 

E 

Williams  et  al  2005 

W05 

w 

E 

I 

Total  Wave 

28 

Total  Current 

5 

Total  W+C 

6 

Total  Equilibrium 

35 

Total  Transient 

4 

Total  Predictors 

28 

Total  Intercomparisons 

24 
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Table  2  Data-sets  used  in  intercomparisons  of  methods  against  data 


Data-set 

Code 

Wave  (W), 
Current  (C) 

Or  both  (W+C) 

Equilibrium  (E) 
or  transient  (T) 

Lab  (L)  or 
field  (F) 

Straub  1939 

S39 

c 

E 

L 

Bagnold  1946 

B46 

W 

E 

L 

Barton  &  Lin  1955 

BL55 

c 

E 

L 

Manohar  1955 

M55 

w 

E 

L 

Inman  1957 

157 

w 

E 

F 

Vanoni  &  Brooks  1957 

VB57 

c 

E 

L 

Laursen  1958 

L58 

c 

E 

L 

Kennedy  1961 

K61 

c 

E 

L 

Yalin  &  Russell  1962 

YR62 

w 

E 

L 

Inman  &  Bowen  1963 

IB  63 

W+C 

E 

L 

Ackers  1964 

A64 

c 

E 

L 

Kennedy  &  Falcon  1965 

KF65 

w 

E 

L 

Guy  et  al  1966 

G66 

c 

E 

L 

Horikawa  &  Watanabe  1967 

HW67 

w 

E 

L 

Vanoni  &  Hwang  1967 

VH67 

c 

E 

L 

Horikawa  &  Watanabe  1968 

HW68 

w 

E 

L 

Carstens  et  al  1969 

C69 

w 

E 

L 

Allen  1971 

A71 

c 

E 

L 

Davies  1971 

D71 

c 

E 

L 

Mogridge  &  Kamphuis  1972 

MK72 

w 

E 

I. 

Pratt  &  Smith  1972 

PS72 

c 

E 

L 

Dingier  1974 

D74 

w 

E 

F 

Banks  &  Collinson  1975 

BC75 

c 

E 

L 

Fok  1975 

F75 

c 

T 

L 

Dingier  &  Inman  1976 

DI76 

w 

E 

F 

Sleath  1976 

S76 

w 

E 

L 

Lofquist  1978 

L78 

w 

E 

L 

Yalin  &  Karahan  1978 

YK78 

w 

E 

L 

Jopling  &  Forbes  1979 

JF79 

c 

E 

L 

Nielsen  1979 

N79 

w 

E 

L 

Alexander  1980 

A80 

c 

E 

L 

DuToit  1980 

DT80 

w 

E 

L 

Miller  &  Komar  1980 

MK80 

w 

E 

L 

Costello  &  Southard  1981 

CS81 

c 

T 

L 

Du  Toil  et  al  1981 

DT81 

w 

E 

L 

Engel  1981 

E81 

c 

E 

L 

Nielsen  1981 

N81 

w 

E 

F 

Hayakawa  et  al  1983 

H83 

w 

E 

L 

Mantz  1983 

M83 

c 

E 

L 

Lambie  1984 

L84 

w 

E 

L 

Shibayama  1984 

S84 

w 

E 

L 

Skafel  &  Krishnappan  1984 

SK84 

w 

E 

L 

Steetzel  1984 

SZ84 

w 

E 

L 

Tanaka  &  Shuto  1984 

TS84 

W+C 

E 

F  +  L 

Vongvisessomjai  1984 

V84 

w 

E 

L 

Sakakiyama  et  al  1985 

S85 

w 

E 

L 

Yalin  1985 

Y85 

c 

E 

L 
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Table  2  Data-sets  used  in  intercomparisons  of  methods  against  data  (continued) 


Data-set 

Code 

Wave  (W), 
Current  (C) 

Or  both  (W+C) 

Equilibrium  (E) 
or  transient  (T) 

Lab  (L)  or 
field  (F) 

Nieuwjaare/  al  1987 

N87 

w+c 

E 

L 

Ribberink  et  al  1987 

RJ387 

w 

E 

L 

Rosengaus  1987 

R87 

w 

E 

L 

van  Rijn  1987 

VR87 

w 

E 

L 

Vongvisessomjai  1987 

V87 

w 

T 

L 

Kos'yan  1988 

K88 

w 

E 

F  +  L 

Mathisen  1989 

M89 

W 

E 

L 

Nelson  Sc  Smith  1989 

NS89 

c 

E 

L 

Lyne  et  al  1990 

L90 

W+C 

T 

F 

Southard  et  al  1990 

S90 

w 

E 

L 

Nielsen  1992 

N92 

w 

E 

L 

Wiberg  &  Nelson  1992 

WN92 

c 

E 

L 

Baas  1993 

B93 

c 

T 

L 

van  Rijn  et  al  1993 

vR93 

W+C 

E 

L 

Willis  et  al  1993 

W93 

w 

E 

L 

Ribberink  &  Al-Salem  1994 

RAS94 

w 

E 

L 

Villaret  Sc  Latteux  1994 

VL94 

w+c 

E 

L 

Van  Rijn  Sc  Havinga  1995 

VRH95 

w+c 

E 

L 

Li  &  Amos  1998 

LA98 

w+c 

E 

F 

Grasmeijer  Sc  Van  Rijn  1999 

GVR99 

w 

E 

L 

Hume  et  al  1999 

H99 

w 

E 

F 

Traykovski  et  al  1999 

T99 

w+c 

E 

F 

Williams  et  al  1999 

W99 

w 

E 

L 

Khelifa  &  Ouellet  2000 

KO00 

w+c 

E 

L 

Williams  et  al  2000 

WOO 

w 

E 

L 

Hanes  et  al  2001 

H01 

w 

E 

F 

O’Donoghue  Sc  Clubb  2001 

ODC01 

w 

E 

L 

Doucette  2002 

D02 

w 

E 

F 

Faraci  &  Foti  2002 

FF02 

w 

T 

L 

Thome  &  Bell  2002 

T02 

w 

E 

L 

Grasmeijer  Sc  Kleinhans  2004 

GK04 

w 

E 

F 

Williams  et  al  2004 

W04 

w 

E 

L 

Doucette  &  O'Donoghue  2005 

DO05 

w 

T 

L 

O'Donoghue  et  al  2005 

OD05 

w 

E 

L 

Smith  Sc  Sleath  2005 

SS05 

w 

T 

L 

Williams  et  al  2005 

W05 

w+c 

E 

F 

Total  Wave 

50 

Total  Current 

22 

Total  W+C 

11 

Total  papers 

83 

Total  Equilibrium 

75 

Total  Transient 

8 

Total  Lab 

71 

Total  Field 

14 
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- 
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CO 
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05 
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- 
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- 

N92 

X 

X 

X 

CO 
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X 

X 

X 

CO 

NS89 

X 

- 

MS86 

X 

- 
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X 

X 

X 

CO 

V84 

X 

X 

X 

CO 

TS84 

X 

- 

CD84 

X 

- 
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X 
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X 

X 

X 

h- 
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X 

X 

X 

X 
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CO 
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CO 
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compar¬ 

ison 

s 

if) 

s 

> 

LO 

CO 

> 

o 

05 

N92 

WN92 

B93 

CO 

O) 

CO 

CD 

z: 

on 

> 

a 

2 

a 

X 

§ 

CO 

5 

CD 

CD 

t- 

KOOO 

HOI 

o 

o 

Q 

O 

FF02 

CO 

o 

ll 

Li- 

3 

CL 

CD 

< 

3 

CD 

1 

in 

o 

O 

Q 

OD05 

W05 

Totals 

a> 

a 

cs 

a 

c 

o 

C/5 

03 

Oh 

a 

© 

a 


T3 


'■© 

a> 


u 

© 


T3 

© 

U 

A< 


ro 

a> 

2 

a 

H 


TR  154 


UNCLASSIFIED 


R  2.0 


Prediction  of  Ripple  Properties  in  Shelf  Seas 

Mark  2  Predictor  for  Time  Evolution  UNCLASSIFIED 


4JhR  Wallingford 


u 

QJ 

CL 

Cl 

cl 

c 

a 


cl 

CL 


O 
u 
u 
a j 


-a 

QJ 


S 

T3 

QJ 


a 

Q 


TT 

a> 

3 

H 


S76 

X 

X 

CM 

DI76 

X 

X 

CM 

F75 

X 

- 

BC75 

X 

- 

D74 

X 

X 

X 

X 

X 

to 

PS72 

X 

- 

MK72 

X 

X 
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X 
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X 

X 

X 

CD 
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- 
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X 

CM 
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Table  6  Prediction  methods  for  ripple  height  and  wavelength,  and  their  functional 
dependence 


Method 

Code 

q/A  and  XI A  are  functions  of: 

Nielsen  (1981) 

N81 

¥ 

Grant  and  Madsen  (1982) 

GM82 

e'w,  d. 

Van  Rijn  (1989) 

VR89 

¥ 

Madsen  (1993) 

Ma93 

(e'„/D?) 

Mogridge  el  al  (1994) 

M94 

Ax 

Wiberg  and  Harris  (1994) 

WH94 

A 

Traykovski  et  al (1999) 

T99 

A 

Faraci  and  Foti  (2002) 

FF02 

¥,  Rew,  A 

Grasmeijer  and  Kleinhans  (2004) 

GK04 

¥ 

Table  7 


References  to  prediction  methods  in  text  books  and  comparison  papers  (see  Table  6 
for  codes) 


Reference  source 

N81 

GM82 

VR89 

Ma93 

M94 

WH94 

T99 

FF02 

GK04 

Text  books 

Sleath  (1984) 

X 

Nielsen  (1992) 

X 

Fredsoe  and 

Deigaard  (1992) 

X 

Van  Rijn  (1993) 

X 

X 

Soulsby  (1997) 

X 

X 

Intercomparisons 

Wiberg  and  Harris 
(1994) 

X 

X 

X 

Traykovski  et  al 
(1999) 

X 

X 

X 

Foti  and  Faraci 
(2003) 

X 

X 

X 

X 

X 

X 

Grasmeijer  and 
Kleinhans  (2004) 

X 

X 

X 

O’Donoghue  et  al 
(2005) 

X 

X 

X 
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Wave  ripple  predictors  with  1000  randomly  distributed  input  parameters ,  plotted 
against  A .  (a)  height  t],  (b)  length  X 
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Figure  2  Wave  ripple  predictors  with  1000  randomly  distributed  input  parameters,  plotted 

against  S7.  (a)  height  //,  (b)  length  A 
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Figure  3  Wave  ripple  predictors  with  1000  randomly  distributed  input  parameters,  plotted 

against  0'^  (a)  height  rj,  (b)  length  X 
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Figure  4 


A 


Non-dimensional  ripple  height  >j/A  from  data-base  plotted  against  (a)  A,  (b)  iP, 
(c)  6'w 
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Figure  5 


A 


Non-dimensional  ripple  wavelength  A/A  from  data-base  plotted  against  (a)  A,  (b)  V7, 
(c)  0'w 
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Figure  6  Test  of  Eq  (4.2,  3)  for  ripple  height  against  data-base 


Lab 
♦  Field 
—  Pred 


Figure  7  Test  of  Eq  (4.1)  for  ripple  wavelength  against  data-base 
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Ripple  heights 


D* 

Figure  8  Comparison  of  new  formula  for  equilibrium  ripple  height  (Eq  5.6a) ,  and  B93 
(modified)  formula  (Eq  5.4),  against  data  from  various  sources 

New  data  is  from  Whitehouse  et  al  (1998),  Damgaard  et  al  (2003),  Lauchlan  (2004)  and  Villaret 
(1994) 


Ripple  lengths 


Figure  9  Comparison  of  new  formula  for  equilibrium  ripple  wavelength  (Eq  5.6b),  and  B93 
(modified)  formula  (Eq  5.5),  against  data  from  various  sources 

New  data  is  from  Whitehouse  et  al  (1998),  Damgaard  et  al  (2003),  Lauchlan  (2004)  and  Villaret 
(1994) 
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D  SO  (/i  m)  at  20°C 


Figure  10  Bedform  existence  plot  by  Van  den  Berg  and  Van  Gelder  (1989),  with  proposed 
wash-out  and  sheet-flow  limits  superimposed.  Additional  new  data  from  B93, 
WMS98  and  DSPW03  is  included 
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(a)  ripple  height  evolution,  d50  =  0.095mm  (b)  ripple  height  evolution,  d50  =  0.238mm 


Figure  11  Comparison  of  Eq  (5.16)  for  ripple  height  evolution  with  the  data  of  B93  for  (a) 
dso  =  0.095mm,  (b)  dso  =  0.238mm;  and  of  Eq  (5.17)  for  ripple  wavelength  evolution 
for  (c)  d$o  =  0.095mm,  (d)  dso  =  0.238mm 
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Figure  12 
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Test  of  wave-ripple  predictor  against  Doucette  and  OfDonoghue  (2005)  experiments 
-  ripple  growth  and  decay 
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Figure  13 
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Test  of  wave-ripple  predictor  against  Doucette  and  OfDonoghue  (2005)  experiments 
-  ripple  growth  and  further  growth 
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Figure  14  Test  of  wave-ripple  predictor  against  three  (non-consecutive)  Doucette  and 
O’Donoghue  (2005)  experiments  -  slow  growth;  growth  from  specified  initial 
height ,  length;  decay  from  specified  height,  length 
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Figure  15  Test  of  current-ripple  predictor  against  Whitehouse  et  al  (1998)  r  ever  sing-flume 
tidal  experiments 
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Figure  16 


Test  of  wave-plus-current  ripple  predictor  driven  by  synthetic  wave  and  current  data 
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Figure  1 7  Test  of  wave-plus-current  ripple  predictor  driven  by  field  data  from  Teignmouth 
( Whitehouse ,  2005) 
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Figure  18  Test  of  wave-plus-current  ripple  predictor  driven  by  Teignmouth  data  with  doubled 
depths  and  doubled  current  speeds 
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Figure  19  Test  of  wave-plus-current  ripple  predictor  driven  by  Teignmouth  data  with  bio¬ 
degradation  switched  ON 
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Appendix  A  Algorithm  for  predicting  the  time- 
evolution  of  ripple  height,  wavelength  and  orientation 
in  sandy  sea-bed  sediments  under  the  influence  of 
waves  and  currents 

A.l  Description  of  prediction  algorithm 

The  purpose  of  the  algorithm  is  to  predict:  time-varying  ripple  heights,  lengths  and 
orientations,  for  any  sandy  sediment,  driven  by  time-series  of:  water  depths,  wave 
heights,  periods  and  directions,  and  current  speeds  and  directions.  It  takes  account  of: 
evolution  (“history”),  threshold-of-motion,  and  wash-out  effects,  and  bio-degradation  of 
ripples.  It  covers  both  wave-generated  and  current-generated  ripples,  switching  between 
them  (in  an  evolving  sense)  depending  on  which  forcing  is  dominant. 

A.2  Inputs 

Steady  inputs 

dso  -  median  grain  diameter  of  sediment  (m) 
ps  =  density  of  sediment  (=  2650  kg  m"3) 
p  =  density  of  water  (kg  m"3) 
v  =  kinematic  viscosity  of  water  (m2  s*1) 
g  =  acceleration  due  to  gravity  (-  9.81  m  s‘2) 

[z  =  height  above  bed  of  current  measurements] 

Swb  =  switch  for  bio-degradation  (=  0  for  OFF,  =  1  for  ON) 

T./Jkb  -  half-life  for  biological  ripple  height  decay  (hours) 

r|b  =  residual  bio-roughness  (m).  (The  average  “trough  to  crest”  height  of 
biologically  induced  bed  features.) 

Time-vatying  inputs  -  time-series  of: 
h  =  water  depth  (m) 

Uw  =  amplitude  of  near-bed  wave-induced  orbital  velocity  (ms-1) 

T  =  wave  period  (s) 

<pw  =  direction  waves  come  from  (degrees  clockwise  of  North) 

U  [or  U(z)]  =  depth-averaged  current  speed  [or  measured  value  at  height  z]  (m  $*') 

(pc  —  direction  currents  go  towards  (degrees  clockwise  of  North) 

Uw  and  T  are  representative  of  regular,  sinusoidal  waves.  O’Donoghue  et  al  (2005) 
showed  from  comparative  laboratory  tests  that  in  regular,  asymmetric  waves  the  best 
representation  for  Uw  is  given  by  Umax,  the  maximum  (usually  onshore-directed)  velocity 
(under  the  wave  crest),  and  for  irregular,  asymmetric  or  symmetric  waves  the  best 
representation  for  Uw  is  given  by  U i/io?  the  mean  of  the  highest  one-tenth  velocities.  In 
irregular  waves,  the  peak-period  Tp  gives  the  best  representation  of  T. 

The  wave  can  alternatively  be  input  as  time-series  of  wave  height  and  period,  where  (for 
the  above  reasons)  the  best  choice  is  Uy\o  and  T i/i0  (mean  height  of  the  highest  1/10  of 
waves,  and  their  corresponding  mean  period).  In  this  case,  the  orbital  velocities  are 
calculated  from  the  depths,  heights  and  periods,  using  linear  wave  theory  or  an 
approximation  to  it. 

The  current  speeds  can  either  be  input  as  depth-averaged  values  (assumed  in  Sections 
A. 4  and  A. 6  below),  or  as  measurements  at  a  fixed  height  z  above  bed.  In  the  latter 
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case,  a  small  modification  to  the  algorithm  is  made  (the  latter  version  is  implemented  in 
the  spreadsheet). 


Initial  values  of  the  ripple  height,  wavelength  and  orientation  are  needed  as  inputs. 
r|0  =  initial  ripple  height  (m) 

Xo  =  initial  ripple  wavelength  (m) 

cpr>o  =  initial  ripple  orientation  (degrees  clockwise  of  North) 


In  most  cases,  the  choice  of  initial  values  will  not  be  important,  after  an  initial  transient 
period.  However,  in  the  case  of  sub-threshold  wave  and  current  conditions  the  initial 
ripples  will  remain  “frozen”,  and  in  conditions  slightly  above  threshold  the  values  will 
change  only  slowly.  Typical  input  values  are  r|0  -  0.02m,  Xq  =  0.2m,  (pf)o  -  shore- 
normal  direction  (wave-dominated)  or  shore  parallel  direction  (current-dominated). 


A.3  Outputs 

Time-series  of: 

r|  =  ripple  height  (m) 

X  =  ripple  wavelength  (m) 

cpr  =  ripple  orientation  (degrees  clockwise  of  North) 

(A.l) 
(A. 2) 


A.4  Calculate  derived  parameters 

Steady  parameters 

s  =  ps/p  =  density  ratio 

g(s-l)lX 


D*  = 


v 


*50 


=  dimensionless  grain-size 


Check  D*  is  in  the  range  of  validity  of  the  algorithm,  1.2  <  D*  <  14,  (approx  d50  in  range 
0.06mm  to  0.7mm).  If  not,  the  algorithm  is  not  applicable. 

z0=d5O/12  (A.3) 

0cr  =0.3/(l  +  l  .2D*)+0.055[1-  exp(-0.02D*)]  -  threshold  Shields  parameter 

(Soulsby  &  Whitehouse  formula)  (A.4) 

Tb  =  3600T./Jfb/ln(2)  -  exponential  time-scale  for  bio-degradation  (s)  (A. 5) 

Time-varying  parameters  (calculated  at  each  time-step) 

If  necessary,  calculate  Uw  from  H  and  T  (or  Hj/jo  and  T 1/10)  using  linear  wave  theory,  or 
an  approximation  to  it.  Then  calculate: 


A  =  UwT/(2ti)  =  amplitude  of  near-bed  wave  excursion  (A. 6) 

A  =  A/d50  (A. 7) 

U 2 

T  =  — -  =  wave  mobility  parameter  (A. 8) 

fw  =  grain-related  wave  friction  factor  (Swart  formula) 

fw  =0.3  for  A  <3.93  (A.9a) 

fw  =  0.00251  exp  [5.21(0.4A)'°19]  for  A  >  3.93  (A.9b) 
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0', 


cD  = 


50 


.  UJJl 

g(s  - 1  y, 

0.40 

ln(h/z0)-l 


=  amplitude  of  oscillatory  skin-friction  Shields  parameter 


-i2 


c  u2 

0/  _ _ ^DU 


x50 


g(s-l)d 

Or  [for  currents  measured  at  height  z]: 


due  to  waves 

=  drag  coefficient  for  depth-averaged  current 
=  skin-friction  Shields  parameter  due  to  current 


c  = 


0.40 


ln(z/z0) 

Q>  _  CzU(z)2 
g(s-l)d 


50 


(A.  10) 
(A.  11) 

(A. 12) 


(A. 13) 
(A.  14) 


A. 5  Calculate  equilibrium  ripple  height,  wavelength  and  orientation  due  to 
wave  forcing,  and  rate-of-change  parameter  p 

Time-varying  parameters  (calculated  at  each  lime-step) 

Check  if  waves  are  above  threshold: 


If  0'w  <  0cr,  ripples  take  pre-existing  values  of  r)  and  X 

If  0'w  >  0cr,  equilibrium  values  of  q  and  X  are  given  by  Eqs  (A.  15),  (A.  16)  and  (A.  17). 
Calculate  equilibrium  quantities: 


^  =  1  + 1 .87  x  1  0’3.a(i  -  expj-  (2.0  x  10~4  a)  5 })]"’ 

A 

(A.  15) 

^3-  =  0. 1 5[l  -  exp{-  (5000/  A)35 }] 

Xeq 

(A.  16) 

Oeq  _  Kq  Tlcq 

(A.  17) 

A  A  ^cq 

(pr,eq  —  *pw 

(A.  18) 

Calculate  rate-of-change  parameter: 

p  =  2.996\|/,O7/(21700  +  y'07) 

(A.  19) 

A.6  Calculate  equilibrium  ripple  height,  wavelength  and  orientation  due  to 
current  forcing,  together  with  Pnand  px. 

Time-varying  parameters  ( calculated  at  each  time-step) 

Check  if  current  is  above  threshold: 


If  0'c  <  0cr,  ripples  take  pre-existing  values  of  q  and  X 
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If  0'c  >  0CT,  equilibrium  values  of  t]  and  X  are  given  by  Eqs  (A. 20)  and  (A.21). 
Calculate: 


Tlmax  =  d5o  •  202  D.° 554  for  1 .2  <  D.  <  1 4 

(A. 20) 

A™*  =  d50  (500  +  1881  D/15)  for  1.2  <D.  <  14 

(A.21) 

0'wo  =1.66D;13  for  D.>  1.58 

(A. 22) 

0'sf  =  2.26  D;1 3  for  D.  >1.58 

0'wo  =  0.916  and  0'sf  =1.25  for  D.  <  1.58 

(A. 23) 

(A. 24) 

Then  calculate  equilibrium  quantities 


qeq  ~  pre-existing  value  for  0  <  0'c  <  0cr 

(A. 25) 

peq  ^lmax  for  0cr  ^00^0  wo 

(A. 25a) 

^leq  =  n  max  ‘  tor  0'wo  <  6'c  <  0'sf 

(A.26) 

Oeq  =  0  fOr6'c>0;f 

(A.  27) 

^eq  =  ^max 

(A.28) 

^Pr,eq  tpc 

(A. 29) 

Calculate  rate-of-change  parameters 


np  _  ^Imax'^max 

‘"[g(s-i)4F 

(A. 30) 

_  20(e;-e„),s 
"  2.5+(e;-e„),s 

(A- 31) 

p  I2(e;-e„),! 

Px-~2.5+(e;-0cr)15 

(A. 32) 

A.7  Set  values  of  coefficients  of  a*  and  bj  for  time-step  tj,  based  on  wave  or 
current  dominance 

Time-varying  parameters  (calculated  at  each  time-step) 

For  ripple  height: 

•  if  0\v  >  0'c,  then: 
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»(0=^-no 


T 


b(t)  =  -£-  +  —  Sw 


T  Tb 

but  if  0'c  >  0'w,  then: 


b(t)=YL  +  ^-Swb 


For  ripple  wavelength: 
•  if  0'w  ^  0'c,  then: 
a(0=|^ 


•  but  if  0'c  >  0'w,  then: 


b(0  =  | 


For  ripple  orientation: 

•  if  0'w  >  0'c,  then: 

b(0  =  Y 

•  but  if  0'c  >  0'w,  then: 


b(0=| 


Set  aj  =  a(tj)  and  b;  =  b(t;)  at  the  i’th  time-step. 


(A. 33) 

(A. 34) 

(A. 35) 

(A. 36) 

(A. 37) 

(A. 38) 

(A. 39) 

(A.40) 

(A.41) 

(A.42) 

(A. 43) 

(A. 44) 
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A. 8  Calculate  Xw  using  Runge-Kutta  integration  (with  over-ride)  for  each 
of  ripple  height,  wavelength  and  orientation 

The  general  equation  to  be  solved  is: 

^-  =  a(t)-b(t)-x(t)  (A. 45) 

dt 


where  x  represents  either  height,  wavelength  or  orientation.  The  time-varying 
coefficients  a(t)  and  b(t)  are  given  in  discrete  form  by  aj  —  a(t,)  and  bj  =  b(t,)  at  the  i’th 
time-step. 

Fourth-order  Runge-Kutta  integration  is  used  to  solve  Eq  (A. 45)  for  each  parameter  in 
turn  (height,  wavelength,  orientation)  by  the  following  steps: 

Calculate  time-step  At  =  t  \+i  -  t  \  (does  not  have  to  be  constant  for  all  t  j).  (A. 46) 


•  If(b.At)  <  1.6,  calculate  in  turn: 

k,  =At.(aj  -b,.xj) 

( 

V 


k2  =  At. 


T(a|  +aI+1)-T(bi  +  bi+1). 


f  1  ,  ^ 

X;  +  — lC 


V 


JJ 


k3  =  At. 


^(ai+a1+i)-^(b,  +  b,+i)- 


r 


x-  +—  ko 
2  ‘ 


\\ 


J) 


k4  =At-(a,+|-b,+i-(x, +  k3)) 


k  =  —  (k,  +  2k2  +  2k3  +  k4) 
6 


x1+,  =  x,+k 


But  if  (b.At)  >  1.6,  instead  calculate: 


xoi  =  xi  + 


vT'x', 


(l  -  exp(-  b,At)) 


(A. 47a) 

(A. 47b) 

(A. 47c) 

(A.47d) 

(A.47e) 

(A.47f) 

(A. 48) 


The  values  of  x  correspond  to  the  predicted  time-series  of  (physical)  ripple  heights  T|, 
wavelengths  X,  and  orientation  cpr. 

Finally,  if  Swb=  1  (bio-degradation  ON),  then: 

height  r|  =  max  (physical  height  r|,  biological  roughness  T|b)-  (A. 49) 

The  predicted  wavelengths  X,  and  orientation  tpr,  are  assumed  to  be  unchanged  by  the 
biological  activity. 
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Appendix  B  Example  input  pages  to  Spreadsheet 


INPUT  PARAMETERS 


|  Enter  values  in  yellow  cells  only 

Computed  values  in  green  cells 

Fixed  values  in  orange  cells 

Name  of  study 
Pre-set  parameters 


Teignmouth  continuous  monitoring  Bursts  15500  -  17000 


Accn  due  to  gravity  = 

9.81 

m/s2 

Sediment  density  = 

2650 

kg/m3 

Physical  parameters 

Water  temperature  = 

18.6 

deg  C 

Water  salinity  = 

32 

ppt 

Water  density  = 

1023 

kg/m  3 

(using  SandCalc) 

Kinematic  viscosity  = 

1.09E-06 

m2/s 

(using  SandCalc) 

d50  (0.06  -  0.7mm)  = 

0.145 

mm 

Height  current  measured  = 

1.21 

m 

Biological  parameters 

Bio-degradation  On  (1);  Off  (0) 

0 

Bio-degradation  half-life  = 

50 

hours 

Residual  bio-roughness  = 

0.005 

m 

For  illustration  purposes  only 

Initial  ripple  properties 

Set  multipliers  to  1  for  normal  use 

Height  = 

0.02 

m 

Multiplier  for  depths  = 

1 

Wavelength  = 

0.2 

m 

Multiplier  for  wave  heights  = 

1 

Orientation  (deg  clockwise  of  N) 

100 

deg  N 

Multiplier  for  currents  = 

1 

D*  = 

g(s-1)d50  = 

[g(s-1)d50A3]A0.5  = 

C_D  for  current  at  height  z 
Threshold  Shields  param.  = 

3.429323279 

0.002262293 

6.89672E-06 

0.001206826 

0.062294672 

(-) 

m2/s2 

m2/s 

(-) 

(-)  (using  SandCalc  Soulsby  &  Whitehouse  method) 

d50(m)  = 

Bio-degradation  timescale  (s)  = 
Bio-degradation  timescale  = 

0.000145 

259685.1074 

72.13475204 

m 

s 

hours 
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